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I.  INTRODUCTION 


The  computer  program  presented  here  calculates  the  Fourier  coef¬ 
ficients  [1,  Eqs.  (87 ) — (89) ]  of  the  electric  and  magnetic  currents  on 
the  object  shown  in  [1  ,  Fig.  1]  excited  by  an  obliquely  incident 
0-polarized  plane  wave.  The  I's  and  V's  in  [1,  Eqs.  (87 ) — (89 ) ]  are 

the  elements  [1,  Eq.  (68a)]  of  X  obtained  by  solving  [1,  Eq.  (65a)] 

n 

for  non-negative  values  of  n  with  B  given  by  [1,  Eq.  (101)]. 

n  " 

The  computer  program  consists  of  a  main  program  and  the  sub¬ 
programs  YZA,  BLOG,  PLANE,  DECOMP,  SOLVE,  and  PRNT.  The  subroutine 
YZA  calculates  the  matrix  elements  [1,  Eqs.  (35)  and  (36)]  which  appear 

in  expression  [1,  Eq.  (33)]  for  the  moment  matrix  T  .  T  is  needed 

n  n 

in  [1,  Eq.  (65a)].  The  subroutine  YZA  calls  the  function  BLOG.  The 

subroutine  PLANE  calculates  V*"®,,  -V^ , ,  -VC  f , ,  and  V^f,  which,  according 

ni  ni  ni  ni 

to  [1,  Eq.  (103)],  determine  of  [1,  Eq.  (101)].  The  subroutine  DECOMP 

n 

decomposes  T  into  the  product  of  a  lower  triangular  matrix  with  an  upper 
n 

triangular  matrix.  The  subroutine  SOLVE  uses  these  triangular  matrices 

to  calculate  the  solution  X  to  [1,  Eq.  (65a)].  The  subroutine  PRNT 

n 

uses  the  elements  of  X  of  [1,  Eq.  (08a) ]  to  calculate  the  Fourier  coef- 

n 

ficients  [1,  Eqs.  (87)— (89) ] .  The  subroutine  PRNT  also  prints  out 
[1,  Eqs.  (87 ) — (89) ] . 

The  main  program  calls  the  subroutines  YZA,  PLANE,  DECOMP,  SOLVE, 
and  PRNT  in  order  to  calculate  and  print  out  [1,  Eqs.  (87 )— (89) ] .  The 
main  program  is  rather  long  because  it  has  to  rearrange  the  Z's  and  Y's 
calculated  by  the  subroutine  YZA  and  the  V's  calculated  by  the  subroutine 
PLANE.  Unfortunately,  the  Z's  and  Y's  come  out  of  the  subroutine  YZA 


i 


r - "i 

2 

arranged  not  as  in  [1,  Eqs  .  (35)  and  (36)]  but  as  on  the  right-hand  sides 

of  [1,  Eqs.  (47)  and  (50)]  for  consecutive  values  of  i'  and  j'.  Likewise,  j 

the  V's  exit  the  subroutine  PLANE  arranged  not  as  on  the  left-hand  sides 
of  [1,  Eq.  (103)]  but  as  on  the  right-hand  sides  of  [1,  Eq.  (103)]  for 

consecutive  values  of  i'.  The  main  program  has  to  rearrange  these  Z's,  j 

Y's,  and  V's  so  as  to  realize  the  transition  from  WrT,  and  JS7,  to  W^T  i 

- ni  nj  - ni 

q± 

and  J  . . 

-nj 

II.  THE  SUBROUTINE  YZA 

The  subroutine  YZA(M1,  M2,  NP,  NPHI,  NT,  IN,  RH,  ZH,  X,  A,  XT, 

AT,  Y,  Z)  calculates  the  elements  of  matrices  Y'  and  Z  defined  by  1 

n  n 
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The  quantities  on  the  right-hand  sides  of  (3)  and  (4)  are  those  on  the 

right-hand  sides  of  [1,  Eqs.  (50)  and  (47)]  with  the  ±  notation  omitted 

and  with  i'  and  j*  replaced  by  i  and  j.  The  <  >  notation  in  (3)  and  (4) 

denotes  the  same  symmetric  product  as  in  [1,  Eqs.  (50)  and  (47)].  The  ± 

notation  in  [1]  serves  to  distinguish  the  surface  (S+  +  A)  from  the 

surface  (S  +  A) .  The  ±  notation  is  not  needed  in  (3)  and  (4)  because  the 

subroutine  YZA  deals  with  only  one  surface  of  revolution  at  a  time.  Later 

on  in  Section  V,  the  main  program  calls  YZA  twice,  once  for  (S+  +  A)  and 

once  for  (S  +  A) .  For  a  particular  value  of  n,  Y^  of  (1)  is  stored  by 

columns  in  Y  and  Z  of  (2)  is  stored  by  columns  in  Z.  Y'  is  placed 
n  n+1  r 

immediately  after  Y^  in  Y.  Similarly,  Z  ^  follows  Z^  in  Z.  Y  and  Z  are 
the  only  output  arguments  of  YZA.  The  rest  of  the  arguments  of  YZA  are 
input  arguments. 

The  input  arguments  of  YZA  have  the  same  meaning  as  those  of  the 

subroutine  YZ  presented  in  [2,  pages  65-79].  In  comparison  with  YZA, 

the  subroutine  YZ  of  [2]  puts  Y  in  Y  and  Z  of  (2)  in  Z  where 

n  n 


Y  = 
n 


!  Ytc 
n 


(5) 


n  n 

The  ijth  elements  of  the  submatrices  on  the  right-hand  side  of  (5)  are 


given  by 


{VS.,  =  -  <Wr.,  n  x  H  ( JS .  ,  0> 
nij  ni  —  -  -nj 


r  =  t,cf) 


s  =  t,< 


(6) 


where  [2,  Eq.  (4)] 


n  = 


Ha 


X 


u 

-t 


(7) 


As  defined  by  (7),  n  is  a  unit  vector  normal  to  the  surface  of  revolution. 


« 

J 

.*a4 


Although  we  have  been  omitting  the  vector  designation  from  vectors  inside 
the  symmetric  product,  we  decided  to  designate  n  and  H  as  vectors  in  (6) 
to  clearly  indicate  that  the  vector  product  n * H  is  intended  there. 

YZ  of  [2]  allows  for  the  Ampere's  law  contribution  to  the  magnetic 
field  H  in  (6).  The  Ampere's  law  contribution  to  the  magnetic  field  is 
the  contribution  due  to  the  value  of  the  electric  current  at  the  field 
point.  The  Ampere's  law  contribution  in  the  subroutine  YZ  of  [2]  is  con¬ 
trolled  by  means  of  the  input  argument  IN  of  YZ  of  [2],  The  action  of  IN 
is  described  in  [2,  page  39]  under  the  assumptions  that 

1)  n  is  given  by  (7) 

2)  n  points  outward  from  the  surface  of  revolution. 

However,  the  action  of  IN  can  be  described  without  recourse  to  assump¬ 
tion  2)  in  the  following  manner. 

IN  =  1,  Magnetic  field  evaluation  on  the  side  of  the  J 

electric  current  sheet  where  the  tail  of  n  is.  / 

IN  =  0,  No  Ampere's  law  contribution.  \  (8) 

IN=-1,  Magnetic  field  evaluation  on  the  side  of  the  ! 

electric  current  sheet  where  the  head  of  n  is.  \ 

In  (8) ,  n  is  visualized  as  piercing  the  sheet  of  electric  current  which 
produces  the  magnetic  field.  Apparently,  assumption  2)  for  n  can  be 
eliminated  from  [2]  by  replacing  the  description  [2,  page  39]  of  the 
action  of  the  input  argument  IN  of  YZ  of  [2]  by  (8). 

Presumably,  the  magnetic  field  H  in  (3)  contains  no  Ampere's  law 

+ 

contribution  because  the  magnetic  field  H  in  [1,  Eq.  (50)]  had  no 
Ampere's  law  contribution.  However,  we  decided  to  allow  for  the  Ampere's 


law  contribution  associated  with  H  in  (3)  in  order  to  enhance  the  use¬ 
fulness  of  YZA.  This  Ampere's  law  contribution  is  controlled  by  means 

of  the  input  argument  IN  of  YZA.  The  action  of  IN  is  described  by  (8) 

/ 

where  it  is  assumed  that  n  is  given  by  (7).  We  do  not  require  n  to 
point  outward  in  (8). 

Most  of  the  statements  in  YZA  are  exactly  the  same  as  those  in  YZ 
of  [2].  In  the  remaining  part  of  this  section,  the  differences  between 
YZA  and  YZ  of  [2]  are  pointed  out  and  explained  and  then  YZA  is  listed. 
Unless  stated  otherwise,  all  line  numbers  and  statement  numbers  cited 
henceforth  in  this  section  refer  to  the  listing  of  YZA. 

In  line  4,  more  space  is  allotted  to  Y  and  Z  than  in  YZ  of  [2]. 
Nevertheless,  the  minimum  allocations  given  in  [2,  pages  65-66]  apply  to 
YZA.  Line  9  in  YZ  of  [2]  has  no  counterpart  in  YZA  because  UG  and  UH 
are  not  used  in  YZA.  On  the  other  hand,  line  17  has  no  counterpart  in 
YZ  of  [2].  Line  17  enables  YZA  to  calculate  the  moment  matrix  for  the 
E-field  solution  [3]  for  the  two  conducting  bodies  in  [1,  Fig.  13].  This 
E-field  solution  was  obtained  by  connecting  the  generating  curves  ABC 
and  DEF  in  [1,  Fig.  13]  to  obtain  the  single  curve  ABCDEF  and  then  de¬ 
leting  all  matrix  elements  associated  with  the  connecting  line  CD.  Line  17 
sets  kp  =  1  at  the  midpoint  of  the  connecting  line  CD.  Here,  k  is  the 
propagation  constant  and  p  is  the  distance  from  the  z  axis.  If  line  17 
were  absent,  this  kp  would  be  zero  and  divisions  by  zero  would  occur. 

The  value  1  assigned  to  kp  in  line  17  is  not  critical  because,  as  stated 
earlier,  all  matrix  elements  associated  with  the  connecting  line  CD  are 
deleted.  Line  17  also  enables  YZA  to  calculate  the  moment  matrix  for 
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an  H-field  solution  for  the  two  conducting  bodies  in  [1,  Fig.  13]. 
However,  this  H-field  solution,  being  a  solution  to  [2,  Eq.  (1)]  without 
the  nx  ,  would  fail  according  to  the  discussion  in  the  third  from  the  last 
paragraph  in  [1,  Section  III].  The  H-field  solution  shown  in  [1, 

Fig.  15]  was  obtained  by  using  the  subroutine  YZ  of  [2]  modified 
by  setting  kp  =  1  at  the  midpoint  of  the  connecting  line  CD  in  [1, 

Fig.  13]. 

The  major  difference  between  YZ  of  [2]  and  YZA  lies  in  the  fact 

that  there  is  an  n*  in  (6)  for  Y^A  calculated  by  YZ  of  [2]  but  not 

in  (3)  for  Y'A  calculated  by  YZA.  Because  VJr .  in  (3)  is  a 
nij  — ni 

tangential  vector,  (3)  can  be  written  as 


Y'^S  *  -  <n  *  Wr.,  n  x  H(JS . ,  0)> 
nij  —  -nx  -  -  -nj 


r  =  t ,  4> 


s  =  t,c 


In  view  of  (7)  and  the  assumptions  that  has  only  a  t  component  and 

vA  has  only  a  component,  (9)  becomes 
— ni 


where 


v,tS  =  - 
nij  nij 


,<Ps  =  'ts 
nij  nij 


s  =  t,<)> 


A®  <vA,  n  x  H(J®.,  0)> 

nij  — ni  — - nj 


(10a) 


(10b) 


r  =  t  ,< 


s  =  t  ,< 


where 


(A  =  u  oA  •  u.) 

— ni  — t  — ni  ^ 

<A  =  (A  ’  0 

— ni  —9  — ni  — t 


(12a) 


(12b) 
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t  S  (f)s 

According  to  (10)-(12),  Y'..  is  the  negative  of  'i  . .  modified  from  (6) 

nij  nij 

and  Y'f®  is  YC®  modified  from  (6).  The  expression  "modified  from  (6)" 
nij  nij 

means  (6)  with  the  t  and  4>  components  of  interchanged. 

Now,  the  matrix  elements  (6)  are  given  by  [2,  Eq.  (12)]  and  it  is 
evident  from  [2,  Eq.  (9)]  that 


t  Ti(t)  -jn<}> 

W  .  •  u  =  -  e 

— ni  — t  p 


(13a) 


4>  Pi(t)  -jn<j) 

W  ,  •  u,  =  -  e 

— ni  — $  p . 


(13b) 


Hence,  (10)  becomes 


where 


..tt  = 

_  0<t>t 

nij 

nij 

= 

YCt 

nij 

nij 

. ,  _ 

nij 

nij 

,»<M>  = 

nij 

nij 

P,(t) 


0f  [2,  Eq.  (12b)]  with  — -  replaced  by  — ^ 


T.(t) 


nij  nij 


(14a) 

(14b) 

(14c) 

(14d) 

(15a) 


of  [2,  Eq.  (12a)]  with 

nij  nij 


T±(t) 


replaced  by 


Pi(t) 


(15b) 


i  i  ..  P.(t)  T.(t) 

Y  j.  =  Y  .  .  of  [2,  Eq.  (12d) ]  with  -  replaced  by  - 

nij  nij  H  P,  F  P 


(15c) 


T  (t)  P  (t) 

YCf .  =  Y  j.  of  [2,  Eq.  (12c)]  with  -  replaced  by  — -  ( 1 5d ) 

nij  nij  x  p  Pi 


Obviously,  the  interchange  of  functions  in  (15)  must  be  accompanied  by 


\  '  SZk  w. 
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appropriate  changes  in  the  limits  of  integration  with  respect  to  t  in 
[2,  Eq.  (12)].  Also,  because  the  magnetic  fields  in  [2,  Eq.  (12)]  as 
it  stands  are  evaluated  on  the  side  of  the  electric  current  sheet  where 
the  tail  of  n  is,  the  Ampere's  law  contributions  to  [2,  Eq.  (12)]  should 
be  multiplied  by  the  input  argument  IN  of  YZA.  These  contributions  are 
the  single  integrals  with  respect  to  t  in  [2,  Eq.  (12)]. 

Consequently,  the  Ampere's  law  contributions  to  (14)  are  given  by 


Y»4>t 

nij 


TT  (IN) 
Pi 


Ci+1 


P.(t)  T  (t)dt 


t . 
l 


-,t4>  _  tt (in) 

nxj  P. 


t 


j+1 

Tj.(t)  P  (t)dt 


t 


j 


(16a) 


(16b) 


where  the  horizontal  bar  on  the  left-hand  sides  of  (16)  denotes  Ampere’s 
law  contribution.  If  the  qth  interval  is  defined  to  be  (tq,  tq+^) ,  then 
the  contributions  to  (16)  due  to  the  integrations  over  the  qth  interval 
are  given  by 


tt  (IN) 
P„ 


P  (t)  T . ( t)dt 

q  j 


j  =  q-i,  q 
j  *  o 

j  +  np  -  1 


niq 


tt(IN) 

P„ 


I  T1(t)  Pq(t)dt 
t 

q 


i  =  q-l,  q 
i  +  0 

i  +  NP  -  1 


(17a) 


(17b) 


where  NP  is  one  of  the  input  arguments  of  YZA.  The  dot  on  the  left-hand 
sides  of  (17)  denotes  Ampere's  law  contribution  due  to  integration  over 
the  qth  interval.  Thanks  to  the  definitions  [2,  Eqs.  (5)  and  (6)],  (17) 


reduces  to 
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A  TT(IN) 

Y'^  =  -9- 

nqj 


2P. 


niq 


A  tt(IN) 

_s _ 


where 

A  =  t“  -  tf 

q  q+i  q 


(18a) 


(18b) 


(19) 


Expressions  (18)  are  the  Ampere's  law  contributions  to  the  matrix  ele- 

r  b 

ments  of  (14)  due  to  the  integrations  over  the  qth  interval. 

The  contributions  to  (14)  referred  to  in  this  paragraph  and  the 
next  paragraph  are  exclusive  of  the  Ampere's  law  contributions.  If  the 
pqth  region  is  the  region  for  which 


t  < 
P  “ 


t  <  t 


P+1 


t  <  t' 

q  - 


< 


t 


q+l 


then  the  contributions  to  (14)  due  to  the  integrations  over  the  pqth  region 

are  given  by  expressions  similar  to  [2,  Eq.  (18)].  If  the  generating  curve 

of  the  surface  of  revolution  is  assumed  to  be  a  series  of  straight  line 

segments  connecting  the  points  t^,  t~,...  tNp  and  if  the  integrations  with 

respect  to  t  are  approximated  by  sampling  the  integrands  at  t  =  t  and 

P 

multiplying  by  A^,  then  the  contributions  to  (14)  due  to  the  integrations 
over  the  pqth  region  reduce  to 
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Since  [2,  Eq.  (22)]  is  equal  to  [2,  Eq.  (24)],  (20)  can  be  re¬ 
written  as 

|  (UC  -  UD)  ,  p  +  0,  q  i  0 
n,p-l,q-l  2 


*  1 1 
yl  LL 

np,q-l 


(UC  -  UD)  ,  p  +  NP-1,  q  )l  0 


*  tt 

Y »  L  L 

n,p-l,q 


-  ^  (UC  +  UD)  ,  p  *  0,  q  #  NP-1 


$,tt 

npq 


j  (UC  +  UD)  ,  p  +  NP-1,  q  +  NP-1 


where 


np,q-l 


*.4>t 

npq 


n,p-l,q 


npq 

$,<M> 

npq 


2  (UA) 
2  (UB) 


\  (UF) 


=  -  f  (UF) 


2  (UE) 


,  q  *  o 
,  q  *  NP-1 
,  P  +  0 
,  p  t  NP-1 


UA  =  right-hand 

side 

of 

[2,  Eq. 

(24a)] 

for 

j  = 

q— i 

UB  =  right-hand 

side 

of 

[2,  Eq. 

(24a)] 

for 

j  = 

q 

UC-UD  =  right-hand 

side 

of 

[2,  Eq. 

(24b)] 

for 

j  = 

q-f 

UC+UD  =  right-hand 

side 

of 

[2,  Eq. 

(24b)] 

for 

j  = 

q 

UE  =  right-hand 

side 

of 

[2,  Eq. 

(24c)] 

UF  =  right-hand 

side 

of 

[2,  Eq. 

(24d) 1 

without 

the 

6  term 

pq 


(24a) 

(24b) 

(24c) 

(24d) 

(24e) 
(24  f ) 
(24g) 
(24h) 
( 24  i) 

(25a) 

(25b) 

(25c) 

(25d) 

(25e) 


( 2  5  f ) 
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The  notation  UA,  UB,...UF  on  the  left-hand  sides  of  (25)  is  the  same 
as  the  notation  used  in  YZ  of  [2].  Expressions  (24)  are  the  contribu- 
tions  to  the  matrix  elements  of  (14)  due  to  the  integrations  over  the 

pqth  region.  No  Ampere's  law  contributions  are  included  in  (24).  The 
Ampere's  law  contributions  are  given  by  (18). 

The  Ampere's  law  contributions  in  (18)  are  different  from  those 
in  [2]  which  consist  of  (2,  Eq.  (23)]  and  the  6  term  in  [2,  Eq.  (24d)]. 
Hence,  the  statements  which  realize  the  Ampere's  law  contributions  in  YZA 
are  different  from  those  in  YZ  of  [2].  Lines  50-51,  74-75,  86-88,  93-99, 
and  424-429  of  YZ  of  [2]  obtain  the  Ampere's  law  contributions.  Line  50 
of  YZA  puts  tt(IN)  of  (18)  in  PN1.  Line  73  puts  the  right-hand  side  of 
(18a)  in  PI. 

The  index  JQ  of  DO  loop  15  obtains  the  subscript  q  which  appears 
in  (18a)  and  (24) .  The  index  IP  of  DO  loop  16  obtains  the  subscript  p 
which  appears  in  (24).  The  purpose  of  the  statement  IF(IP.NE.JQ)  in 
lines  121,  270,  and  424  of  YZ  of  [2]  is  to  check  if  the  IPth  interval 
coincides  with  the  JQth  interval.  If  the  IPth  and  JQth  intervals  coincide, 
then  the  numerical  evaluation  of  the  integrals  G  in  [2,  Section  III]  is 
affected  and  Ampere's  law  contributions  are  taken  into  account.  However, 

YZ  of  [2]  was  not  designed  to  accommodate  a  generating  curve  which  closes 
upon  itself. 

The  subroutine  YZA  allows  for  a  generating  curve  which  closes  upon 
itself.  Such  a  curve  is  treated  in  the  following  manner.  A  triangle 
function  whose  peak  is  at  the  first  data »point  is  needed  and  is  obtained 
by  overlapping  the  last  interval  of  the  generating  curve  with  the  first 
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interval  as  in  [1,  Fig.  5],  This  triangle  function  is  called 

In  addition  to  TNp_2(t),  the  process  of  overlapping  obtains  the  pulse 

function  P.tT)  ,(t).  P._n  ,  (t)  is  not  wanted  because  it  is  identical  to 
NP-I  NP-1 

P^(t).  The  effect  of  P^Tp  ^(t)  can  be  eliminated  by  discarding  all  matrix 

r  s 

elements  (3)  and  (4)  for  which  either  W  .  or  J  ,  contain  P„_  ,  (t).  These 

~ni  — nj  NP-1 

elements  must  be  discarded  after  exit  from  YZA,  because  YZA  does  not  con¬ 
tain  any  logic  for  discarding  matn  elements. 

Unfortunately,  (18)  was  derived  for  a  generating  curve  which  does 
not  overlap  on  itself.  If  the  last  interval  of  the  generating  curve 
overlaps  the  first  interval,  the  quantities  on  the  right-hand  side  of 
(18)  are  still  correct  but  more  values  of  i  and  j  are  needed  in  (18) 
when  q  *  1  and  q  =  NP-1  to  allow  for  the  overlapping  so  as  to  account  for 
all  Ampere's  law  contributions  to  the  matrix  elements  (3).  It  will  be 
shown  that  all  Ampere's  law  contributions  can  be  accounted  for  by  de¬ 
fining  ZIP  to  be  the  electrical  distance  from  the  center  of  the  JQth  interval 
to  the  center  of  the  IPth  interval  and  using  the  statement  IF(ZIP.NE.O.) 
instead  of  the  statement  IF(IP ,NE. JQ) .  The  statement  IF(ZIP.NE.O.)  appears 
in  lines  110,  259,  and  408.  If  the  generating  curve  does  not  overlap  on 
itself,  the  action  of  the  statement  IF(ZIP.NE.O. )  is  the  same  as  the  action 
of  the  statement  IF(IP.NE.JQ)  and  all  Ampere's  law  contributions  are  accounted 
for.  However,  if  the  last  interval  of  the  generating  curve  overlaps  the 
first  interval,  then  it  remains  to  be  shown  that  (ZIP  =0.)  obtains  all 
Ampere's  law  contributions. 

If  the  last  interval  of  the  generating  curve  overlaps  the  first 
interval,  then  ZIP  is  zero  not  only  for 
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IP  -  JQ  =  m. 


m=l,2, . . .NP-1 


Case  1 


but  also  for 


Case  2 


and  for 


IP  =  NP-1 
JQ  =  1 


J 


Case  3 


r  s 

In  what  is  to  follow,  the  matrix  elements  Y*  ,  of  (3)  are  viewed  as 

nij 

interactions  between  testing  functions  and  expansion  functions.  The  combi¬ 
nation  of  Cases  1,  2,  and  3  will  give  all  Ampere's  law  contributions  if  it 
covers  all  possible  interactions  between  parts  of  testing  functions  on  the 
mth  interval  and  parts  of  expansion  functions  on  the  mth  interval  for 
m  =  1,2, . . .NP-2. 

With  regard  to  the  Ampere's  law  contributions  (18),  the  last  state¬ 
ment  is  more  general  than  necessary  because  not  all  interactions  between 
testing  functions  and  expansion  functions  are  involved  in  (18).  For 
example,  the  interactions  Y^j  between  t  directed  testing  functions  and  t 
directed  expansion  functions  are  not  involved  in  (18)  and  neither  are  the 
interactions  Y^^  between  (j)  directed  testing  functions  and  <f>  directed 

expansion  functions.  However,  the  matrix  elements  .  and  Y^f .  of  (6) 

nij  nij 

calculated  by  YZ  of  [2]  do  have  Ampere's  law  contributions.  It  would  be 
useful  to  know  that  YZ  of  [2]  can  be  modified  to  allow  for  a  closed  gen¬ 
erating  curve  by  using  the  same  technique  as  in  YZA.  This  technique  con¬ 
sists  of  overlapping  the  last  interval  of  the  generating  curve  with  the 
first,  replacing  the  statement  IF(IP.NE.JQ)  by  the  statement  IF(ZIP.NE.O.) , 


and  deleting  the  matrix  elements  associated  with  PNp  ^(t)  after  exit 
from  the  subroutine. 

It  is  evident  that  the  portion  of  Case  1  for  which 
IP  =  JQ  =  m,  m  =  2....NP-2 


covers  all  possible  interactions  between  parts  of  testing  functions  on 
the  mth  interval  and  expansion  functions  on  the  rath  interval  for 
m  *  2,3,...NP-2.  However,  clarification  is  needed  with  regard  to  the 
first  interval  because,  due  to  the  overlapping,  the  first  interval  is 
sometimes  disguised  as  the  (NP-l)th  interval. 

All  possible  interactions  between  parts  of  testing  functions  on 
the  first  interval  and  parts  of  expansion  functions  on  the  first  interval 
are  listed  as 


D  TNP-2(t)  Wlth  TNP-2(t) 

2)  T*(t)  with  T^p_2(t) 

3)  P1(t)  with  T^p_2(t) 
A)  T^p_2(t)  with  T*(t) 

5)  T*(t)  with  T*(t) 

6)  Px(t)  with  T*(t) 

7)  TNP-2(t)  Wlth  Pl(t) 

8)  lj(t)  with  Px(t) 

9)  Pp( t)  with  P1(t) 


In  each  of  the  foregoing  9  interactions,  the  first  function  is  the  part 
of  the  testing  function  and  the  second  function  is  the  part  of  the  expan¬ 
sion  function.  t!L  0(t)  is  the  second  part  of  the  (NP-2)th  triangle 
NP-2 
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function,  the  downward  sloping  part.  T^(t)  is  the  first  part  of  the  first 
triangle  function,  the  upward  sloping  part. 

The  interactions  covered  by  Cases  2  and  3  and  the  parts  of  Case  1  for 
which  IP  =  JQ  =  1  and  IP  =  JQ  =  NP-1  are  listed  as 


2) 

T*(t) 

with 

TNP-2(t) 

3) 

PL(t) 

with 

TNP-2(t) 

10) 

T*(t) 

with 

PNP-l(t) 

ID 

PL(t) 

with 

PNP-l(t) 

IP  =  1 
JQ  =  NP-1 


Case  2 


4)  TNP-2(t)  With  Tl(t) 


12)  PNP-l(t)  Wlth  Tl(t) 


7)  Tjjp_2(t)  with  Px(t) 


I»  =  NP-1 
JQ  =  1 


Case  3 


13)  pNP_1<t)  with  p1<t) 


5)  T*(t)  with  T*(t) 


6)  P1(t)  with  Tx(t) 


8)  T*(t)  with  P1(t) 


9)  Px(t)  with  Px(t) 


IP  =  JQ  =  1 


Part  of  Case  1 


U  TNP-2(t)  With  TNP-2 ' 

14)  PNP-l(t)  Wlth  TNP-2 * 

15)  TNP-2(t)  Wlth  PNP-1* 

16)  PNP-1^  Wlth  PNP-1* 


IP  =  JQ  -  NP-1 


Part  of  Case  1 
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▼ 


The  preceding  16  interactions  were  numbered  so  as  to  facilitate  compari¬ 
son  with  the  9  interactions  in  the  last  paragraph. 

Of  the  16  interactions  in  the  last  paragraph,  numbers  1  to  9  are  the 
interactions  in  the  second  from  the  last  paragraph  and  numbers  10  to  16 
are  to  be  discarded  because  they  contain  P^p  ^(t).  Drawing  an  equivalence 
between  the  interactions  in  the  last  paragraph  and  the  interactions  in  the 
second  from  the  last  paragraph,  we  conclude  that  the  combination  of  Cases  1, 
2.  and  3  covers  all  possible  interactions  between  parts  of  testing  functions 
on  the  mth  interval  and  parts  of  expansion  functions  on  the  mth  interval 
for  m  =  1,2,...  NP-2.  Hence,  (ZIP  =  0.)  obtains  all  Ampere's  law  contribu¬ 
tions  when  the  last  interval  on  the  generating  curve  overlaps  the  first 
interval.  If  such  overlapping  is  intended,  then  the  input  arguments  RH  and 
ZH  of  YZA  must  satisfy 


RH(NP-l)  =  RH(1) 


RH(NP) 


RH(2) 


ZH(NP-l)  =  ZH(1) 


ZH(NP)  =  ZH(2) 


(26) 


exactly.  Otherwise,  the  computed  values  of  ZIP  will  not  be  zero  in  Cases 
2  and  3. 

With  the  intention  of  showing  that  YZA  implements  (24),  we  list  in 
Table  1  some  variables  in  YZA  whose  values  are  different  from  those  in 


YZ  of  [2]. 


Table  1.  Comparison  of  some  variables  in  YZA  with  those  in 
YZ  of  [2]. 


Line 
in  YZA 

Variable 
in  YZA 

Expression  in  terms  of 
variables  in  YZ  of  [2] 

379 

W1 

2 .  *W1 

380 

W2 

2.*W2 

385 

H1C 

-.5*H1C 

387 

H3C 

2.*H3C 

388 

H2C 

2 . *H2C 

389 

H3C 

2 .  *H3C 

390 

W3 

- . 5*W3 

391 

W4 

- . 5*W4 

392 

W5 

-.5*W5 

405 

UC 

2.*UC 

406 

UB 

2.*UB 

407 

UF 

- . 5*UF 

411 

UA 

2 .  *UA 

412 

UB 

2.*UB 

415 

UC 

- . 5* (UC-UD) 

416 

UD 

- . 5* (UC+UD) 

417 

UE 

2 .  *UE 

In  Table  1,  the  variable  in  the  second  column  is  defined  by  the  statement 
in  YZA  whose  line  number  is  given  in  the  first  column.  The  expression  in 
the  third  column  is  the  value  of  the  variable  in  the  second  column.  The 
values  of  the  variables  appearing  in  the  third  column  are  the  values  which 
exist  in  YZ  just  after  the  counterpart  to  the  variable  in  the  second  column 
has  been  defined  in  YZ.  For  example,  the  value  of  W1  in  the  first  line  of 
the  third  column  is  defined  in  line  391  of  YZ  of  [2],  As  an  exception,  the 


19 


value  of  UD  intended  in  the  third  from  the  last  line  of  the  third  column 
is  that  defined  by  line  433  of  YZ  of  [2].  Table  1  was  constructed  under 
the  assumption  that  no  Ampere's  law  contributions  come  into  play.  This 
amounts  to  assuming  that  lines  425-429  of  YZ  of  [2]  are  not  executed  and 
that  lines  409-410  of  YZA  are  not  executed.  The  Ampere's  law  contribu¬ 
tions  (18)  will  be  considered  later. 

The  last  6  variables  in  the  second  column  of  Table  1  are  used  to 

construct  the  matrix  elements  Y',®  of  (3).  The  statements  in  YZA  which 

nij 

perform  this  task  are  listed  in  Table  2. 

rs 

Table  2.  Construction  of  matrix  elements  Y' 

nij 


Line  Statement  Matrix 

in  YZA  in  YZA  element 


461,476 

Y(K1) 

= 

Y(K1) 

+ 

uc 

466,483 

Y(K2) 

= 

Y(K2) 

+ 

uc 

448,477 

Y(K3) 

= 

Y(K3) 

+ 

UD 

453,484 

Y(K4) 

= 

UD 

458,471 

Y(K5) 

= 

Y(K5) 

+ 

UA 

445,472 

Y(K6) 

= 

UB 

449,462,478 

Y(K7) 

= 

Y(K7) 

+ 

UF 

454,467,485 

Y(K8) 

= 

UF 

489 

Y(K9) 

cs 

UE 

Y’tt 

n,p-l,q-l 

Y.tt 

n,p,q-l 

Y,tt 

n,p-l,q 

Y.tt 

npq 

Y  »4>t 
np,q-l 

Y  »4>t 
npq 

Y,t<J> 

n,p-l,q 

Yit<J> 

npq 

Y»<W 

npq 


In  Table  2,  the  variable  being  modified  by  the  statement  in  the  second 
column  represents  the  matrix  element  in  the  third  column.  For  example, 
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Y(K1)  which  appears  in  the  first  line  of  the  second  column  of  Table  2  is 

the  storage  location  of  the  matrix  element  Y,tt:  ,  , .  The  subscripts  n.p, 

n, p-1 , q-1  K 

and  q  appearing  in  the  third  column  of  Table  2  are  given  in  terms  of  vari¬ 
ables  in  YZA  by 

n  =  M  +  Ml-1 
P  =  IP 
q  =  JQ 

Here,  M  is  the  index  of  DO  loop  31,  Ml  is  one  of  the  input  arguments  of 
YZA,  IP  is  the  index  of  DO  loop  16,  and  JQ  is  the  index  of  DO  loop  15. 

In  view  of  the  fact  that  the  variables  UA,  UB,  UC,  UD,  UE,  and  UF  appear¬ 
ing  in  the  second  column  of  Table  2  are  given  by  the  last  6  entries  of 
the  third  column  of  Table  1,  it  is  now  evident  that  YZA  implements  (24). 

From  (18),  the  correct  values  of  the  Ampere’s  law  contributions  are 

i(A^Tr(IN))/(2p^) .  In  the  paragraph  which  follows  (25)  we  established  that 
line  73  stores  (A^tt(IN))/ (2p^)  in  PI.  In  the  paragraph  prior  to  the  intro¬ 
duction  of  Table  1,  we  concluded  that  the  Ampere's  law  contributions  come 
into  play  whenever  ZIP  =  0.  Accordingly,  lines  409-410  are  executed  when 
ZIP  =  0.  The  effect  of  line  409  is  to  add  PI  to  the  values  of  UA  and  UB 

calculated  by  lines  411  and  412.  Line  410  subtracts  PI  from  the  value  of  UF 

calculated  by  line  407.  The  previously  mentioned  variables  UA,  UB,  and  UF 
appear  in  the  second  column  of  Table  1  and  were  transferred  to  the  second 
column  of  Table  2.  It  is  now  evident  from  Table  2  that  execution  of 

lines  409-410  adds  PI  to  Y'^  ,  and  Y’^)t  and  subtracts  PI  from  Y,t:^  , 

np,q-l  npq  n,p-l,q 

and  Y,t:^.  Hence,  YZ  takes  proper  account  of  the  Ampere's  law  contributions, 
npq 
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001 C  LISTING  OF  THE  SUBROUTINE  YZA 

002  C  THE  SUBROUTINE  VZA  CALLS  THE  FUNCTION  BLOG 

003  SUBROUTINE  YZAi Ml . M2. NP. NPHI . NT . I N »RH. Z H. X. A . XT  * AT, V . Z > 

0  04  COMPLEX  Y I  2209  I . Z 1 2209 » .U ,U I .U2.U3.U4.U5.U6 .HI A . H2 A. H3 A, GA C  48 ) 

005  COMPLEX  G8(48).GC( 48) .GDI  48) .GE< 48 » .HIB .H2B .H3B.HIC .H2C.H3C.H4A 

006  COMPLEX  H5A.H6A.H4B.H58. H6B.UA, ua.UC.UO.UE.UF.GIAI 10).G2A< I 0> 

00  7  COMPLEX  G3  A(IO),GtB(tO) . G2  B ( I  0 ) . GJO  { I0).G1C(10) .G2CI  10) ,G3C(  I0> 

0  08  COMPLEX  G4A( 101.G5AI 1 Ol .G6A( I  01 ,G4B( 10) .G58( I  0 > *G6B(  I0>  .CMPLX 

009  OIMENSION  RH 143 1 . ZHI 4  3  I ■ X ( 4  8 ) , A< 4 8 ) . XT{  10).AT(  10). RSI  42)  .ZSI42) 

0  10  OIMENSION  0142 ) . OR I  42 ) . OZ 1 42  )<DM(42  )«Ct (48 ) . C2 ( 48 ) ,C3( 48) .C4(  200) 

Ott  OIMENSION  C5I200) . C 6 < 200) >R2( 10) .Z2( I O) *R7( 1 0) . Z7 (  10) 

012  CT~2. 

013  CP=.l 

014  OO  10  1*2 «NP 

015  12*1-1 

016  RS(I2>*.5*IRHf D+RHII2)) 

017  IF(RS( (2I.E0.0.)  RS(12)=l. 

018  ZS (I2>=»5*( ZHI I ) 4ZH< 12)) 

019  Ol  *.5*<  RH(  I  >— RH<  12  )  ) 

020  02*.S*(ZH( 1 )-ZH€ 12) ) 

021  D( I2)=Sa«T(01 *D1*D2*02) 

022  0RII2)=01 

023  DZ(I2J*02 

024  OM( 1 2 )*DI t  2) /RSI  12) 

025  10  CONTINUE 

026  M3 » M2 -Ml  41 

027  M  4— M  1-  1 

028  PI  2*1 .570796 

029  PP*9. 869604 

030  OO  1 1  K* 1 . NPHI 

031  PH=PI24(X(K)4l .) 

032  C1(K)*PH 

033  C2IK)*PH4PH 

034  SN*SINI .5*PHI 

035  C3(K)*4.4SN4SN 

036  Al*PI24AIK> 

037  D4*.S 4 A 14C3IK ) 

038  DS*A14COSIPH) 

0J9  06=A I 4SINI PH I 

940  MS*K 

041  OO  29  M*l.M3 

042  PHM*(M4+M«*PH 

043  A2*COS(PHM) 

044  C4(MS)*04*A2 

C4S  C5(MS)*OSAA2 

046  C6IM5 )=06*SIN<PHM) 

047  M5*M5+NPHI 

048  29  CONTINUE 

049  II  CONTINUE 
050  PN1*3« 1415934IN 

051  U*10..t.) 

052  U1*.S«U 

053  U2*2.*U 

054  MP=NP— 1 

055  MT*MP-1 

056  N*MT4MP 

057  N2 N*NT4N 

058  N2*N«N 

CS9  JM*-1-N 

060  00  IS  JO* I • MP 
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061 

KQ=2 

062 

IF(JQ.EQ.l)  K0=l 

063 

IF< JQ.EQ.MP)  K0=3 

064 

R l=RS( JO) 

065 

Z1=ZS(3Q) 

066 

01=0(30) 

06/ 

D2=OR( JQI 

068 

03=02(30) 

069 

04=02/Rl 

070 

05=01 /Rt 

071 

S V=D2/D 1 

072 

CV=D3/D1 

073 

Pl=PNl*OS 

0  74 

P3=2 - *D 1 

0  75 

P 4=2. *D4 

0  76 

P5=04*D4 

077 

P 6=01*01 

078 

P7=P6*01 

079 

T6=CT*D 1 

080 

T62=T6+Ot 

08  1 

T62=T62*T62 

082 

R6=CP*R 1 

083 

R62=R6*R6 

oe* 

DO  12  1_  =  1.NT 

085 

06=XT  ( 1_  ) 

086 

R2(L)=R1 *02*06 

087 

Z2CL  J=Z1+D3*06 

088 

12  CONTINUE 

0  89 

OO  16  IP=1.MP 

090 

R3=RS(  IP) 

091 

Z3=2S< IP) 

092 

R4=Rl— R3 

093 

24=21-23 

094 

U3=02*U1 

095 

U4=03*U1 

096 

DO  40  L=l.NT 

097 

07=>R2(L  )— R3 

098 

08=22  <  L )— 23 

099 

R7(l.)=R3*R2<l.l 

too 

Z  7(L  )=07*07*08*08 

10  1 

40 

CONTINUE 

1  02 

PH=R4*SV+Z4*CV 

103 

A l=ABSC PH) 

104 

A2=  AOS  <  R4  *CV-Z4*SV I 

105 

06=  A2 

1  06 

IF(Al.EE.Dl)  GO  TO  26 

107 

06=41-01 

1  08 

06= SORT ( 06*06  *A2*A2 > 

l  09 

26 

21  P=R4*R4  ♦•24*24 

1  10 

IF (ZIP.NE.0..AND.1R6. 

1  1  1 

Z5=2IP 

1  12 

R5=R3*R| 

1  13 

PHM=.5*R3*SV 

1  14 

OO  33  K«t.NPHI 

1  IS 

A 1=C 3( K ) 

1  16 

RR«Z5+RS*AI 

117 

HI A=0. 

1  18 

H2A=0. 

119 

H3A=0» 

120 

H4A=0. 

121 

122 

123 

124 

125 

126 

127 

128 
I  29 

130 

131 

132 
1  33 
I  34 
1  35 
136 
1  J7 
1  38 
1  39 

140 

141 
1  42 
143 
1  44 
145 
I  46 

147 

148 

149 
1  50 

151 

152 

153 
1  54 

155 

156 

157 

ise 

159 

160 
161 
162 
163 
l  64 

165 

166 
167 
166 
169 
1  70 
171 
1  72 
173 
1  74 
I  75 
176 
1  77 
1  78 
1  79 
1  80 


HSA=0. 

IMRR.LT.r62)  GO  TO  34 
DO  35  L=1.NY 
M>Z7(U8)7<U«*I 
R=SQRT<W) 

SN=-SINIR) 

CS=COS(R) 

06— AT IL l/R 

HI 8=D6/M*CMPLX(CS-R*SN,SN*R*CS> 

HIA-H1B »H 1  A 
H2 8—  XT  (  L  )  •  H  1  8 
H2A=H2B*H2 A 
H3A=XT< L) 4  H284H3A 
H4B=  064CNPL  X(CS.SN) 

H4A=H4B*H4 A 
HSA=XTI L) *H4B*H5A 
35  CONTINUE 
GO  TO  36 
34  OO  37  L= 1 .NT 

*=Z7IL  )  8)7(11*41 
R-SQRT ( W) 

IFIR.GT..5)  GO  TO  14 

CS=R*(M*(.6944444E-2-»*. 1 736 1 1IE-3I-. I 25) 

SN=9* I . 3333333 E— l— 9*  » 1 190476E-2 )-. 3333333 
H1B=AT(L) *CMPL  X(CS# SN ) 

CS=R*lir*(  .  4  1 66  66  7E  —  l  —  .  1  388  889E-2*  M  )  -  .  5  ) 

SN=»MW*(  .  1984 126E-3**-. 83333 JJE-2 ) ♦. 1666667)- 1. 
H48=AT(L)*CMPLX(CS.SN) 

GO  TO  43 
14  SN— — 5  I NIR) 

CS=COS(R) 

06=AT(L)/R 

HI  8=064 <  (CMPLXICS-R*SN,SN*R*CS)- l.)/M-.5) 
H4B=06«CMPLX ICS-1 • *  SN ) 

43  HIA=H1B*M1A 
H28=XT<L)*HI8 
H2A=H2B*H2A 
H3A=XT( L) *H2B*H3A 
H4A=H4B*M4A 
H5A=XT(L) *H4B*H5A 

37  CONTINUE 
A1=PH*PHM*A1 
A2  =  A0SI A  1 ) 

RsRR-A2*A2 
06=42-01 
07=42*01 
062=06* 06 
072=07*07 

O  8=  SORT  1 062*RI 
09= SORT  I O  72*R ) 

IF  I R— I RR* 1 .E-5 ) )  52.52,53 

52  IFI06.LT. 0.)  STOP 
K4=.5/D62— .5/072 
GO  TO  54 

53  K*=<07/'09-06/0al/R 

54  IFIO6.GE.0.)  GO  TO  38 
M=ALOG( 107*09) *<-C6*08l/R) 

GO  TO  39 

38  9=ALOGI  107*09)  /I  06*03)  ) 

39  lll  =  tW4*.S*V)/Dl 
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isi 
i  a  2 

I  dJ 
I  84 

185 

186 
18/ 
1  88 

189 

190 

191 

192 

193 

194 

195 

1  96 

197 

198 

199 

200 
201 

2  02 

203 

204 

205 

206 
2  07 
208 

209 

210 
211 
212 

213 

214 

215 

216 

217 

218 
2  19 
220 
22  I 
222 

223 

224 

225 

226 

227 

228 

229 

230 

231 

232 
2  33 

234 

235 

236 

237 

238 

239 

240 


85=A  2/D  I 

82=1.  5*1  09—  08  )  -  1  ./D9  ♦  l  ./D 8 1/P6-M5 *8  l 

83  =  1.254  <0  7*09— D6 *08 ) »8-R*( *4 ..25*8 I l/P7-85*(2. *82.8S*8 1 ) 
84=8/01 

8S=<0»-08-A2*8  »/P6 
1F1A1.GE.0. I  GO  TO  27 
82=-82 
8S=-8S 

27  HI A=8t  *H1 A 
H2A=82»H2A 
H3A=83*H3A 
H4A=84*H4A 
H5A=85»H5A 
36  GAIK l=H! A 
G8(K)=H2A 
GC(K>=H3A 
GO(KI=H4A 
GE1K1=HSA 
33  CONTINUE 
K  1=0 

OO  45  tt=l.M3 
H1A=0. 

H2A=0. 

H3A=0. 

HI  8=0. 

M2B=0. 

H38=0. 

H|C=0. 

H2C=0. 

H3C=0. 

H4A=0. 

HSA=0. 

H6 A=0. 

H48=0. 

H58=0. 

H6B=0. 

OO  46  K= 1 • NPMI 

K1=K1*1 

06=C41KI) 

07=C51K1» 

D8=C6(  K 1 ) 

U  A=G  A  1 K  > 

U0=G8(K) 

UC*GC«K) 

UO=GO<K) 

UE=GE(K) 

H1A=06*UA*H|A 
N2A=07*UA*H2A 
N3A=OS*UA*H3A 
H 10=06*  U8*H 18 
H2B«07*UB.H2B 
H3B=08*US*M3B 
H1C=06*UC*HIC 
H2C*074UC*H2C 
H3C=08*UC»H3C 
H4A«06*U0*H4A 
H5A»074U0*H5A 
H6 A=08*U0*H6A 
H48=064Ue*H48 
H5B=07*ue*H58 


♦•v  ' 
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241  H6a*D6» UE4-H68 

242  46  CONTINUE 

243  G1A(MI=H1A 

244  G2A(MI=H2A 

245  G3A(M)-H3A 

246  G18(MI=HIS 

247  G2S(M)=H2B 

248  G3B( M  1=H3B 

249  GIC(M)=mC 

250  G2C«M)-H2C 

251  G3CtM»=M3C 

252  G4A(M)=H4A 

253  GSA(M)3H5A 

254  G6A(M*=H6A 

255  G4S<M»=H4B 

256  GSB(M|-H58 

257  G6Q<  M  1=H68 

258  45  CONTINUE 

259  IMZIP.NE.O.I  GO  TO  47 

260  41=05405 

261  08*0. 

262  09=0. 

263  DO  63  K=I.NPHI 

264  D8=084A(KI/SQHT(C2<K) »At I 

265  09=094- A  (  K 1 4BLOG  (OS/Cl  (  K  )  ) 

266  63  CONTINUE 

267  42=3.141593/05 

268  O8=IBCOGtA2l-PI240«»/<«54RI } 

265  09= 2. /R I  41 BLOG< A2I 4A2*8L0G(  I ./A2) )- 3. 1 4  1 593/0 1409 

270  DO  67  M=1.M3 

271  G1A(MI=084G1 A(M> 

272  G2A(M)=0. 

273  G2B<MJ=0. 

274  G2C(M)=0. 

275  G3A(MI*0. 

276  G5A<MI=094GSA(MI 

277  67  CONTINUE 

278  GO  TO  47 

279  41  OQ  25  M=t«M3 

280  GIA(Mt=0. 

281  G2AtN)=0. 

282  G3A{ N  >  =  0. 

283  Gie(N)=0. 

284  G2B(M)=0. 

285  G38<M|=0. 

236  GIC<M»=0. 

287  G2CCN»=0. 

238  G3C(M»=0. 

289  G4A<N»=0. 

290  G5A<MI=0. 

291  G6A(M)=0. 

292  G4B(M)=0. 

293  G5B ( M 1=0. 

294  G6B<M»*0. 

295  25  CONTINUE 

296  00  13  L=l • NT 

297  RS=R7  <C> 

298  25=271 L ) 

299  OO  17  K=I.NPHI 

300  4=254RS4C3IK» 


sot 

302 

303 

304 

305 
3  06 

307 

308 
300 

310 

311 

312 

313 

314 
3  1 S 

316 

317 

318 
310 

320 

321 

322 
3  23 

324 

325 

326 

327 

328 

329 

330 

331 

332 

333 

334 

335 

336 

337 

338 

339 

340 

341 

342 

343 

344 

345 

346 
34  7 
14  8 

349 

350 

351 

352 

353 

354 

355 

356 

357 

358 

359 

360 


R=SORT(tf) 

SN=-SIN(R) 

CS=COS(R> 

GAIK  )=CMPLX<CS-R*SN»SN*R*CS)/<M*R) 

GOCK)=CNPLX(CS.SN)/P 
17  CONTINUE 

IFIR62.LE.ZS)  GO  TO  51 
06=0. 

07-0. 

09^0. 

OO  62  K* 1 • NPHK 
M2=C2<K) 

«=1./IZ5*R5*W2) 

M|*A<K) *SORT<«> 

06=06 *W 1**2 *N 

07=07* W 1*(.5*«*< 1 . 4.1254 W*R5*M2*W2) ) 

09=09 *« I 
62  CONTINUE 
W1=RS/ZS 
M2=PP*Ml 
«*SQRT(W2) 

M3=1.*W2 
R=SORT< M3) 

W4=SQRT (RS) 

•5=AL0G<M*R) 

08=-PI2«06-(W/R-NS> /<RS*U4) 

06* .5 *08 

07=1 1 «/R* ( M-t . 1 25*. 1666667**2 )/»3 ) *. 1 25*M5  >/R5*.S*W5)/W4-P  1 2 *07 
09=MS/«4-P 1 2*09 
51  A1  =  AT(L) 

A2=XT<L)*A1 
A3*XT(L )*A2 
K  1=0 

OO  30  M=1 .M3 
M=M*M4 
HI  A=0. 

H2A*0. 

H3A=0. 

H4A=0. 

HSA=0. 

H6A*0. 

DO  32  K=1 • NPHI 
K1=K1 *1 
HI B=GA<  K) 

W4=C4(KI) 

W5=C5(K1  ) 

M6=C6(KI ) 

HlA=M4*HIO*HIA 
H2A=«S*HIB*H2A 
H3A=M6*H1 B*H3 A 
HIB*60(K> 

H4A=»4*H1B*H4A 
HS A=MS*H1 B*HSA 
H6A=M6*  H l B  *H6A 
32  CONTINUE 

IFIR62.LE.Z5)  GO  TO  44 
H 1 A=06*H1 A 

H2A=D7-<W«M*I.)*D6*H2A 

H3A=V*08*H3A 

H5A*09*HSA 


36  1 

362 

363 

364 

365 

366 

367 

368 

369 

370 

371 
3  72 

373 

374 

375 

376 

377 
3  7e 

379 

380 

381 
3  9  2 
383 

3  84 

385 

386 
397 

388 

389 

390 

39  1 

392 

393 

394 

395 

396 

397 

398 

399 

400 

40  1 

4  32 

403 

404 

405 

406 

407 
40e 
409 

41  0 
41  1 

412 

413 
4  1  4 
415 
4  I  6 
4  1  7 
4  18 

419 

420 


44  G1A(M)=A1*H1 A 40 l AIM) 

G2AIM)=A1*H2A+G2A(M) 

G3AC  N)iAI  4H3A4-G3AI  M  ) 

Gt 8<N)=A2*H1 A*G18(M) 

G2BCM)~A2*H2A+G2B(M) 

G3S( N )=A24H3A+G38(  Ml 
GIC(M)=A3*H1A*G1C(M) 

G2CIM)=A3*H2A+G2CIM» 

G3C< M )=A3*H3A*G3CI M> 

G4A(M)=A1*H4A*G4A(M) 

G5AI M )= At *H5A+G5AI Mt 
G6A( M >=A 14H6A4G6A ( M) 

G4  B ( M 1= A24H4A4G4B  CM) 

GSB( M )  =  A2  *HS A+G5B ( M ) 

G6B( M1=A2*H6A4G6B( M) 

30  CONTINUE 

13  CONTINUE 

47  A2=OIIP) 

M1=A2*(R4*03-Z4*02) 

*2=>-A2*H3*03 

A 3=02 < IP) 

06= DR IIP) 

07=24406 

09=03*06 

H1C= I  02*1 R3*A3*07 )-R 1 *09 ) *U I 

08=A2*0I 

H3C=D8*U2 

H2C=Z4*H3C 

H3C=03*H3C 

*3=0 I » I 07-R4* A3 ) 

*4=01 *< 09— 02*A3) 

*5=— D1 *R 1 *A3 
Al=OR( IP) 

US=AI*U3 
U6=A3*U4 
06=-02*A2 
O  7=0  I  •  A  I 
A3=0Mf IP) 

JM=  JN 

DO  31  M=1.M3 
H2Ap=  G2AIM) 

HI A=Gt AIM) 

H2B=G2BIM» 

H1B=G1BIN) 

UOb*I*H2A4«2*H1A 

UB=MI*H28+M2*HIB 

UF=M3*<  H2A»04*H2B ) 4*4* I H2B4044 G2C ( M ) ) 4*5* (HI A4P4 *M1 B4P5 *G 1 C ( M ) ) 
IFf 2IP.NE.0. )  GO  TO  48 
UC=UCt-PI 
UF=UF-P I 

48  UA=UC-UB 
U8=UC+UB 
H3A=G3A(M) 

H3B=G38(M) 

UC=H1 C* (H3 A— H3B) 

UO=HlC*< H3A4H3B) 

UE=H2C*f  H3A*D4*H3B >  4HJC  *( H3B 404*G 3C  C M) ) 

H5A=G5AIM) 

HSBaGSB(M) 

H4A*C4A(M) *H5A 
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421 

H4Ss046(Ml 4H5B 

422 

H6A=G6A4M) 

423 

H6B*G6a<M) 

424 

H3A=US«H5A»U6*H4A 

4  25 

H18«U5*H5B«U6*H4a 

426 

HI AxHJA-HIB 

427 

H2A*H3A4HIB 

42e 

H3A=-U14H4A 

4  29 

HI e^D6«H6A 

430 

M*MFM4 

431 

At»««A3 

432 

M2  B- D6  •  H6  9—  A 1  *H4A 

4J3 

H30=D7*4H6A+O4*H68i 

434 

H4A=M*05*H4A 

435 

K1*!P»JN 

436 

K2=K14 | 

437 

K  3=K 1 4N 

438 

K4=X.2  +  N 

439 

K5=K2*MT 

440 

K6=K4+MT 

44  1 

K7=K3*N2N 

442 

K8*K4+N2N 

443 

K9=K8VMT 

444 

GO  70  1 18.20. 191. KO 

445 

18  V(K61=ua 

446 

Z4K6)=H1B«-H28 

447 

IFI  IP.EQ.l )  GO  TO  21 

44e 

Y4K3)=Y4K3J»UO 

445 

VI K7)=V4K7)4UF 

450 

Z1K3)=Z(K3) 4H2A— H3A 

451 

Z ( K71 -l 4 K7 ) 4H3B-H4 A 

452 

IF (IP.EQ.MP)  GO  TO  22 

453 

21  VI K4 >— UD 

454 

Y4K8I=UF 

455 

Z«  K4»=H2A4H3A 

456 

Z1K8I=H38 FH4A 

457 

GO  TO  22 

458 

19  YIK5)=Y(KSJ*UA 

455 

Z4K5>=Z4K5»*H18-H2B 

460 

IF(IP.EQ.I)  GO  TO  23 

461 

V4K1)=V4K1 I4UC 

462 

V4K71=Y4K7»4UF 

463 

ZIKl )zZ(KI )4H1 A+H3A 

464 

Z4K71— Z4K7I 4M38—H4  A 

481 

465 

IF41P.EO.HPi  GO  TO  22 

482 

466 

23  V4K2)*V4K2 1 4UC 

483 

467 

V4K8)*UF 

484 

468 

Z4K2)=Z4K2)*H1 A-H3A 

485 

465 

ZC  *C8  >=H38*H4A 

486 

4  70 

GO  TO  22 

487 

471 

20  V 1 K5 1  =  Y  C  K5 ) 4UA 

488 

472 

V4K6)*U8 

489 

473 

ZIK3I*Z4KS)»HIB-H2a 

490 

4  74 

Z 4  K6 ) sHI B4H2B 

491 

475 

IF4IP.E0.I)  GO  TO  24 

4  92 

4  76 

74X11=74*1 )4UC 

493 

477 

Y ( K3 ) »V  4  K3 ) 4U0 

494 

4  78 

VIK7)*V4K7i 4UF 

495 

479 

Z4K1  >«Z4K  1 14HIA4H3A 

496 

480 

Z  4  K3) *Z  4  K3 ) 4M2  A— H3  A 

497 

Z  4  K  7 1«Z { K7 I ♦H3B-H4 A 
IF4 IP.EQ.MP)  CO  TO  22 
2*  y(K2)*V< K2  I 4UC 
V(K4|sUO 
V4K8i»UF 

Z4K2)sZ4K2)»Ht A-H3A 
Z4K4)»H2A+H3A 
Z4K8>*=H38tH4  A 
22  Y(K91*UE 

Z4K9I«U2*4D8*< H5A4044H5B )-A I *M4AI 
JH* JM*N2 
31  CONTINUE 
16  CONTINUE 
JN*JN»N 
15  CONTINUE 
RETURN 
END 
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III.  THE  FUNCTION  BLOG  AND  THE  SUBROUTINES  PLANE,  DECOMP,  AND  SOLVE 

The  function  BLOG  is  exactly  the  same  as  in  [3,  page  56].  The 
only  difference  between  the  subroutine  PLANE  (of  the  present  report) 
and  the  subroutine  PLANE  of  [3,  pages  57-62]  is  that  the  statement 

IF(Rl.EQ.O.)  R1  =  .5 

in  line  39  of  PLANE  has  no  counterpart  in  the  subroutine  PLANE  of  [3, 
pages  57-62].  The  action  of  this  statement  is  similar  to  that  of  the 
statement  in  line  17  of  the  subroutine  YZA.  Except  for  a  difference 
in  the  actual  space  allocated  to  the  variables  UL,  SCL,  and  IPS,  the 
subroutine  DECOMP  is  the  same  as  the  subroutine  DECOMP  of  [3,  pages  63-64]. 
Minimum  allocations  in  DECOMP  are  the  same  as  those  in  the  subroutine 
DECOMP  of  [3,  pages  63-64].  The  subroutine  SOLVE  differs  from  the  sub¬ 
routine  SOLVE  of  [3,  pages  63-64]  only  in  the  actual  space  allocated  to 
the  variables  UL,  B,  X,  and  IPS.  Minimum  allocations  in  SOLVE  are  the 
same  as  those  in  the  subroutine  SOLVE  of  [3,  pages  63-64]. 


001  c 
002 
003 
004 
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007 
008 
006 
OIC 
Oil 
012 
C 1  3 
0  14 
0  IS 

cie 

017 
018 
019 
0  20 
021 
022 
023 
024 
025 
026 
027 
0  28 
029 
030 
031 
032 
033 
034 
0  35 
036 
037 
033 
C39 
040 
041 
042 
043 
044 
04  5 
046 
047 
048 
049 
0  50 
051 
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055 
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058 
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LISTING  OF  THE  FUNCTION  BLOG 
FUNCTION  BLCG(X) 

IFIX.GT..I)  GO  TO  1 
X2*X4X 

BLOG-*  <  I  .0754X2-.  1666667)  4X24  1.  )  *X 
RE  I URN 

I  SLOG*ALOG( X4SQRTI 1.4X4X)) 

RETURN 

ENO 

LISTING  OF  THE  SUBROUTINE  PLANE 

SUBRUUf INE  PLANEIM1.M2.NF.NP.NT.RH.2H.XT.AT.THR.R) 

COMPLEX  RI240) .U.Ut .UA.UB.FAI 10). FBI  1 0 ) . F2 A. F2B. F 1 A. FI B .U2 . U3. U4 
COMPLEX  US.CMPLX 

DIMENSION  RH(43).ZH( 43>.XT<  10 ) . A T I  1 0) . THR I  3) ,CS<  3) ,SN<3) .R2II  0) 

O I MENS  I ON  22 (10) *0  J ISO ) 

MP*NP— I 

MT*MP— l 

NSMT4MP 

N2*26N 

OO  II  K=  l  .  NF 

X*THR<K) 

CSIKI*COS<XI 
SN<KJ=»SINC  X) 

11  CONTINUE 
U*I0..l.) 

U|3]. 14 15934U44M1 

M3*Ml4t 

M4=*M2  43 

IFCMI.EO.O)  M3»2 

M5=MI42 

M6*M2 42 

OO  12  I  P=  1 .  MP 

K2*IP 

I— IP4 | 

DR*.S4|RHf I )-RH(IP) ) 

02s*.  54 1  ZHI I  )-2H<  IP)  t 
01*S0RT IOR4DR40240Z ) 

R1 3.254 (RH< I I4RHI IP) ) 

IFIRI.EO.O.)  R1-.5 
Zt*.54(2H<  D42HI  IP)  ) 

OR3.S4DR 
D2*DR/ft 1 
OO  13  L*1.NT 
R2 (L)*R1 40R*XT (L ) 

Z2  <L) *2 1 40Z4XT  CL) 

13  CONTINUE 

OO  14  K=1.NF 
CC*CS(K) 

SS*SN(K) 

03*  DR ACC 
04>-0Z4 SS 
05*0 I 4CC 
OO  23  M*M3«M4 
F  A I M) *0 • 

F81M)*0. 

23  CONTINUE 

00  IS  L*1.NT 
X*SS4R2|L) 

IF! X.GT..5E-7)  GO  TO  19 
OO  20  M*MJ.M4 


061 
C62 
063 
064 
065 
066 
067 
066 
069 
0  70 
071 
0  7? 

0  7  J 
0  74 
075 
076 
0/7 
078 
0  79 
080 
081 
082 
083 
084 
C35 
086 
087 
C88 
039 
090 
091 
092 
093 
094 
095 
096 
097 
C98 
099 
1  00 
101 
102 
1  03 
104 
l  05 
I  06 
107 
I  08 
I  09 
t  10 
1  1  t 
I  12 
113 
»  14 
115 

ue 

117 

t  \  a 

t  19 

I  20 


H3(M)a0. 

20  CONTINUE 
SJ(2I»I< 

S»l. 

GO  TO  18 
|9 

IFIX.LT. .5)  M* I  I .8 4AL0G 1 0 ( X I 

IFIM.LT.M4j  M=M4 

8J<MI=0. 

BJC JM)s|. 

OO  16  J=»4.M 
J2aJM 
JMaJM-1 
Jt  a JM— I 

BJ(  JM)=J1/X*SJI  J2)-OJ<  JM*2» 

16  CONTINUE 
S=0. 

IFv'  M.LE.4 1  GO  TO  2* 

OO  17  J=4.M.2 
S=*S»BJ<  J> 

17  CONTINUE 

24  SaBJ I2I42.6S 
Id  Aft  G*22  4  L 1 4CC 

UAaAT(LJ/S*CMPLXICOS(AftG» .SIN( AHG>» 

UBaXTIL )*UA 
OO  25  M=M3.M4 

FAIM»aaJlM*4UA4FAIM» 

FBIM»=BJIM J4UB4FQIM* 

25  CONTINUE 
15  CONTINUE 

IFIMI.NE-OI  GO  TO  26 
FACI 1=-FA! 3* 

F  81 1 ) *—FBC  3» 

26  UAaUI 

OO  27  MaM5.M6 

M7*M-l 

MflaMFI 

F2 A=  UA* IF A ( M8 ) 4F A I M7 ) ) 

F2B*UA4<F81  M0I4F8IM7II 
U8aU«UA 

FI A*UB»<FA<M8»-FAtM7)> 

F 1 8* UB • (FBI  MB  1— FB ( M7>  I 
U«a04*UA 

U2aOJ4FI A4U44FAIM) 

U3=03*F  IB4U44FB1M) 

UAaQft 4F2A 
U5*OH4F2 B 
K 1«K2- I 
K4aK|4N 
X5=K24N 

ft! K24  MT) =— 054 ( F2A402*F2BI 
«OC5»MT  1-0I*(FIA»0  24'IB> 

ifiip.eo.i »  60  T0  ** 
HIKI»*HIKI»4U2-U3 
RIK4»aft(K4 I4U4-U5 
tFUP.EO.Mft)  GO  TO  22 
21  «(X2)*U24U3 
RIX5)*U44U3 
2?  K2*K2*N2 
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121 

ua»us 

122 

27 

CONTINUE 

123 

14 

CONTINUE 

124 

12 

CONTINUE 

125 

RETURN 

126 

ENO 

127C 

LISTING  OF  THE  SUBROUTINE  OECOMP 

128 

SUBROUTINE  OECOMP ( N . IPS .UL I 

129 

COMPLEX  UL16241 ) .PIVOT .EM 

130 

OIMENSION  SCLI 791 .IPSI79) 

131 

OO  S  1*1.  N 

t  32 

(PSUisl 

133 

RN*0. 

1  34 

J1*I 

136 

OO  2  3=1 .N 

136 

ULM=A8S1REAL<UL<  J1 ) ) )  4-ABSi  At  MAGI  ULl  Jill) 

137 

J1=J1*N 

138 

IFIRN-ULM)  1.2.2 

139 

1 

RNsULM 

140 

2 

CONTINUE 

141 

SCLI I 1*1 ./RN 

142 

5 

CONTINUE 

143 

NMl*N— l 

144 

K2=0 

145 

OO  17  K*t.NMt 

146 

BIG=0. 

147 

OO  11  I*K.N 

148 

IPbIPSII) 

149 

IPK=IP+K2 

150 

SIZE*! ABS<REAL<UL< IPKi ) )*ABS< AIMAGIULI IPX) > ) )*SCL< IP) 

151 

IFISIZE-BIG)  11.11.10 

152 

10 

BIG*S1ZE 

153 

<PV*I 

15* 

1 1 

CONTINUE 

181 

IPsIPSI  1) 

1  55 

IFIIPV-K)  14.1S.14 

182 

X(l)*B(IP) 

1  56 

1  4 

J*1PS<K) 

183 

DO  2  1*2. N 

157 

IPS IK  I* IPS 1 1PVI 

184 

IP=IPSI I) 

158 

IPSIIPV)=J 

185 

IP8=IP 

159 

IS 

KPP*IPS1K|4K2 

186 

IMI*I-1 

160 

PIVOT=UL<KPP) 

187 

SUM»0. 

161 

W1*K*1 

188 

OO  1  3*1, IMI 

162 

OO  16  I*KP 1 *N 

189 

SUM=SUM+UL C I P) «X I J ) 

163 

KP*KPP 

190  I 

IP*(P»N 

164 

1P*(PS< 1 ) 4K2 

191  2 

X( I >*8( IPS)- SUM 

165 

EM* -ULI  IP  1  /P  1 VOT 

192 

K2=N*«N-| ) 

166 

i  a 

UL  (  IP)*— EM 

193 

IP*IPS<N)*K2 

1  67 

OO  16  J=KPl.N 

194 

X<N)*X<N)/UL< IP1 

168 

IPslP.N 

195 

OO  4  IBACK*2.N 

169 

KP*KP4N 

196 

1= NP l—l BACK 

170 

ULI  IP»=>ULI  1P)4EM«UL(KP) 

197 

K2=K2-N 

1  71 

16 

CONTINUE 

198 

1P1= IPSI I ) 4K2 

1  72 

K2*K24N 

199 

IP1  *  (  ♦! 

I  73 

17 

CONTINUE 

200 

sum*  o* 

1  74 

RETURN 

201 

ip»ipi 

1  75 

ENO 

20  2 

OO  3  J*IP1 .N 

1  76 

LISTING  OF  THE  SUBROUTINE  SOLVE 

203 

IP-1P4N 

»  TY 

SUBROUTINE  SOLVE ( N . IPS .UL . 8. X > 

204  3 

SUM=SUM«UL 1  IP ) 4X ( J ) 

17« 

COMPLEX  UL <6241 ) .01 79» .XI 79 > . SUM 

205  4 

X(I)*(X(I )— SUM) /ULI IPI ) 

1  79 

OIMENSION  IPS! 791 

zoe 

RETURN 

1  80 

NP 1 *N* 1 

207 

ENO 

% 
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IV.  THE  SUBROUTINE  PRNT 

The  subroutine  PRNT  uses  knowledge  of  the  superscripted  I  's  or 

nj 

V^'s  and  the  k  p"'s  in  a  portion  of  fl,  Eqs.  (87a),  (88a),  and  (89a)] 
and  in  the  corresponding  portion  of  [1,  Eqs.  (87b),  (88b),  and  (89b)] 
in  order  to  calculate  and  print  out  these  portions  of  [1,  Eqs.  (87)- 
(89)].  By  calling  the  subroutine  PRNT  repeatedly,  it  is  possible  to 
print  out  all  the  quantities  on  the  right-hand  sides  of  [1,  Eqs.  (87)- 
(89)]  except  those  for  the  cases  in  which  j  =  l  and  j  =  N+  in  [1,  Eq.  (87a)] 
and  j  =  M  +1  and  j  =  M  +M  in  [1,  Eq.  (89a)].  PRNT  is  not  designed  to 
treat  these  cases.  They  are  included  in  [1,  Eqs.  (87a)  and  (89a)]  merely 
for  convenience. 

Except  for  the  cases  in  which  j*l  and  j  =  N+  in  [1,  Eq.  (87a)] 
and  j  = M+ + 1  and  j  =  M+ +  M  in  [1,  Eq.  (89a)],  any  portion  of  [1,  Eqs.  (87a), 
(88a),  and  (89a)]  can  be  written  as 

Cl  *  XX(J  +  JX)/RA(J  +  J3),  J  «  1,2, . . . J1  (27) 

and  the  corresponding  portion  of  [1,  Eqs.  (87b),  (88b),  and  (89b)]  can 
be  written  as 

C2  *  XX  (J+J1+JX) / (RA(J+J4-1)  +  RA(J+J4)),  J  =  1,2,...J2  (28) 

In  (27)  and  (28),  XX  represents  the  I^'s  or  Vnj's  H-l »  an^  repre¬ 
sents  the  k+  p-’s  of  [1].  In  (27),  Cl  is  either  e  or  2 j .  In  (28), 

n 

C2  is  either  4j  or  2e  .  The  sum  of  the  two  RA's  in  the  denominator  of 

n 

(28)  is  due  to  [1,  Eq.  (12)].  The  complex  numbers  (27)  are  t  directed 

\ 

currents,  either  electric  or  magnetic.  The  complex  numbers  (28)  are  4-’ 
directed  currents.  The  set  of  complex  numbers  (28)  is  linked  to  the  set 
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(27)  by  the  fact  that  the  XX' s  used  in  (28)  occur  immediately  after  those 
used  in  (27). 

First,  the  subroutine  PRNT  (Jl,  J2,  J3,  J4,  RA)  prints  out  the  real 
part,  the  imaginary  part,  and  the  magnitude  of  each  of  the  complex  numbers 
(27)  under  the  heading  which  appears  in  statement  10.  Then,  PRNT 
prints  out  the  real  part,  the  imaginary  part,  and  the  magnitude  of  each  of 
the  complex  numbers  (28)  under  the  heading  which  appears  in  statement  13. 
The  variables  in  (27)  and  (28)  enter  PRNT  by  means  of  the  arguments  of 
PRNT  and  the  statement 

COMMON  Cl,  C2,  XX,  JX 

in  line  6.  If  Jl  =  0,  nothing  from  (27)  is  printed  out.  If  J2  =  0, 

nothing  from  (28)  is  printed  out.  All  the  arguments  of  PRNT  are  input 

arguments.  The  common  variables  Cl,  C2,  and  XX  are  input  variables. 
However,  the  common  variable  JX  functions  as  both  an  input  variable  and 
an  output  variable.  PRNT  adds  J1  +  J2  to  the  original  value  of  JX. 

Minimum  allocations  in  PRNT  are  given  by 
COMPLEX  XX(J2  +  Jl  +  JX) 

DIMENSION  RA(Max (Jl  +  J3 ,  J2  +  J4)) 

where  Max  denotes  the  larger  of  the  two  values  in  the  parentheses.  Of 
course,  the  space  allocated  to  XX  in  the  calling  program  must  be  exactly 
the  same  as  the  space  allocated  to  XX  in  line  4  of  PRNT. 

DO  loop  11  prints  out  the  complex  numbers  (27).  Line  14  obtains 

(27).  DO  loop  14  prints  out  the  complex  numbers  (28).  Line  27  obtains 


(28). 


00 1C  LISTING  OF  THE  SUBROUTINE  PRNT 

002  SUBROUTINE  PRNTI J I . J2. JJ. J4 iRAI 

003  COMPLEX  C1.C2.U 

004  COMPLEX  XXI 791 

305  OIMENSION  RAI43I 

006  COMMON  C1.C2.XX.JX 

007  IFIJI.EO.OI  GO  TO  IS 

006  UR  I TE (3.101 

009  10  FORMAT  1*0  REAL  JT  IMAG  JT  MAG  JT» > 

010  K3 - J3 

Oil  OO  11  J=*1.J1 

012  X3>K341 

013  JX^JX+l 

014  U=C1/RA(K3)*XX( JX) 

015  W=CA8S(U) 

016  WRITEI3.12)  U.U 

Cl  7  12  FORMAT! IX. 3EII.4) 

018  11  CONTINUE 

019  IS  IFIJ2.E0.0)  RETURN 

020  WRITE (3.13) 

021  13  FORMAT! *0  REAL  JP  IMAG  JP  MAG  JP*> 

022  K4  — J4 

023  OO  14  J31.J2 

024  JX=JX+1 

025  X3=K4 

026  K4=K4+1 

027  UxC2/(RA(K3)4RA(K4))*XX( JX) 

028  W=CABS(UI 

029  WRITEI3.12)  U.U 

030  14  CONTINUE 

031  RETURN 

032  ENO 
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V.  THE  MAIN  PROGRAM 

The  main  program  accepts  input  data  and  calls  the  subroutines  YZA, 

PLANE,  DECOMP,  SOLVE,  and  PRNT  in  order  to  calculate  and  print  out  [1, 

Eqs.  (87)  —  (89) ] .  The  input  data  are  read  from  punched  cards  according  to 

READ(1 , 15)  NT,  NPHI 

15  FORMAT(2I3) 

READ(1,10) (XT(K) ,  K=l,  NT) 

READ (1,10) (AT (K) ,  K=l,  NT) 

10  FORMAT (5E14. 7) 

READ(1,10) (X(K) ,  K=l,  NPHI) 

READ (1 ,10) (A(K) ,  K=l,  NPHI) 

READ(1 , 16)  NA,  NB,  MA,  MB,  MC,  LA,  LB,  LC, 

Ml,  M2,  BK,  UR,  ER,  THR(l) 

16  F0RMAT(10I3/4E14 .7) 

READ(1,18) (RA(I) ,  1=1,  NA) 

READ(1 , 18) (ZA(I) ,  1=1,  NA) 

18  FORMAT (10F8. 4) 

READ(1 , 18) (RB(I) ,  1=1,  NB) 

READ(1,18) (ZB(I) ,  1=1,  NB) 

Most  of  the  input  variables  in  the  main  program  represent  variables 
in  [1]  and  [2].  Table  3  relates  input  variables  in  the  main  program  to 
variables  in  [1]  and  [2].  The  input  variable  NT  in  Table  3  represents  both  n^  and 
n^.  It  Is  assumed  that  nfc  =  nT  in  the  main  program.  The  Gaussian  quadrature 

<»t>  <nt)  cy  (n  ) 

data  x?  ,  ,  A^,  ,  x^  r  ,  and  A^  T 


are  given  in  [4,  Appendix  A],  It  is 
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Table  3. 

Input  data  for  the  main 

program. 

Variable  in 
main  program 

Variable  in 
[1]  or  [2] 

Equation  or  figure 
number  in  [1]  or  [2] 

NT 

ni 

t  T 

[2,  Eqs.  (35)  and  (79)] 

NPHI 

[2,  Eq.  (36)] 

XT 

(nt) 

V  1 

,  £’  =  1,2, . . .nt 

[2,  Eq.  (35)] 

AT 

A<v 

V  ' 

,  £'  =  1,2, ..  .n 

[2,  Eq.  (35)] 

X 

<y 

,  £  =1,2,.. .n^ 

[2,  Eq.  (37)] 

A 

(n.  ) 
A  ♦ 

,  £  =  1,2, . . .n^ 

[2,  Eq.  (36)] 

NA 

N+ 

[1,  Fig.  4] 

NB 

N  “ 

[1,  Fig.  5] 

MA 

M+ 

[1,  Fig.  4] 

MB 

M~ 

[1,  Fig.  53 

MC 

M 

[1,  Figs.  4  and  5] 

BK 

k+ 

[1,  Eqs.  (73)— (84 ) ] 

UR 

u~/y+ 

[1,  Fig.  1] 

ER 

e~/e+ 

[1,  Fig.  1] 

THR(l) 

9t 

[1,  Eq.  (92)3 

RA 

-+ 

pr 

j  =  1,2, ...N+ 

[1,  Eq.  (11)3 

ZA 

-+ 

Y 

j  =  1,2, ...N+ 

[1,  Eq.  (11)3 

RB 

V 

j  =  1,2, ...N" 

[1,  Eq.  (11)3 

ZB 

Y 

j  =  1,2, ...if 

[1,  Eq.  (11)3 
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If  LB  =  1,  it  is  assumed  that  M  >.  1  and  that  N  -  M  -  M  >_  2.  LC=1 
implies  that  the  object  under  consideration  reduces  to  a  homogeneous 
material  obstacle  with  parameters  (p  ,  e  )  in  contrast  with  the  parameters 
of  the  external  environment.  As  for  Ml  and  M2,  [1,  Eqs.  (87) -(89)]  are 
printed  out  for  n  =  Ml,  Ml  +  1....M2. 

Minimum  allocations  in  the  main  program  are  given  by 

COMPLEX  XX(N) ,  YP (M*JA*JA) ,  ZP(M*JA*JA), 

YM(M*JB*JB) ,  ZM(M*JB*JB),  R(2*M*JA), 

T (N*N) ,  Y (N) 

DIMENSION  XT (NT) ,  AT (NT) ,  X(NPHI),  A(NPHI),  RA(NA) , 

ZA(NA) ,  RB(NB) ,  ZB(NB) ,  IPS(N) 

where 


4 

4 

N  = 

l 

L1A(I)  +  l  L1B(I)  +  L2(l)  +  L2 (2)  +  L3(l)  +  L3(2) 

(29) 

1=1 

1=1 

M  = 

M2  - 

-  Ml  +  1 

(30) 

JA  = 

2  * 

NA  -  3 

(31) 

JB  = 

2  * 

NB  -  3 

(32) 

The  L's  on  the  right-hand  side  of  (29)  are  calculated  in  the  main  program. 

At  any  rate,  N  is  the  order  of  the  moment  matrix  T^  given  by  f 1 ,  Eq.  (33)]. 

At  the  end  of  this  section,  the  main  program  is  listed  along  with 
sample  input  and  output  data.  The  sample  data  are  for  the  object  of 
[1,  Fig.  13]  with  k+a  =  2.5  and  ^  =  2  as  in  [1,  Fig.  17].  The  output  data 
printed  under  the  heading  "Electric  current  on  first  part  of  outside  con¬ 
ductor"  are,  with  reference  to  [1,  Eq.  (87)], 

j  =  2,3, ...M+  (33) 

j  =  1,2,3, ...M+  (34) 


n  Ujj’ 

4>+0  + 
n  (Zy' 


\ 


AO 


The  real  and  imaginary  parts  and  magnitudes  of  (33)  are  printed  under 
the  heading  "Real  JT  Imag  JT  Mag  JT."  The  real  and  imaginary  parts 
and  magnitudes  of  (34)  are  printed  under  the  heading  "Real  JP  Imag  JP 
Mag  JP. "  The  output  data  printed  under  the  heading  "Electric  current 
on  second  part  of  outside  conductor"  are 


Jt+6(t+),  j  =M++M+1,  M+  +  M  +  2, . . .N+  -  1 
n  j  . 


J^+0(tt),  j  =  M+  +  M,  M+  +  M  +  1 , . . .N+  -  1  -  LA  ( 

n  J 

The  output  data  printed  under  the  heading  "Electric  current  on  first  part 
of  inside  conductor"are,  with  reference  to  [1,  Eq.  (88)], 

J  *  2'3-"-M  1 


r  (tj),  J  -  1.2,3.. ..H  <: 

The  output  data  printed  under  the  heading  "Electric  current  on  second  part 
of  inside  conductor"  are 


Jn  j  =  M  +  M  +  1,  M~  +  M  +  2....N”  -  1 


9 (tT)  »  j  =  M_  +  M,  M-  +  M  +  1,...N~  -  1  -  LB 
n  j 

The  output  data  printed  under  the  heading  "Electric  current  in 
aperture"  are,  with  reference  to  [1,  Eq.  (87)], 


Jt+9(tt),  .1  =  M+  +  1,  M+  +  2,...M+  +  M 

n  J 


J^+9(tt),  j  *  M+  +  1,  M+  +  2,...  M+  +  M  -  1  -  LC 
n  j 


I 
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If,  as  in  (1,  Fig.  18],  there  is  no  conducting  surface  below  the  aper¬ 
ture,  then  j  =  M+  +  1  is  deleted  in  (41).  If,  as  in  [1,  Fig.  10],  there 
is  no  conducting  surface  above  the  aperture,  then  j  =  M+  +  M  is  deleted 
in  (41).  The  output  data  printed  under  the  heading  "Magnetic  current  in 
aperture"  are,  with  reference  to  [1,  Eq.  (89)], 


and 


M^6(tJ),  j  =  M+  +  2,  M+  +  3,...M+  +  M  -  1 


M^0(tj),  j  =  M+  +  1,  M+  +  2,...M+  +  M  -  1  -  LC 


(43) 


(44) 


If,  in  any  of  (33)-(44),  the  upper  limit  on  j  is  less  than  the  starting 

value  of  j,  then  no  values  of  j  are  to  be  taken.  For  example,  if  M+  *  1 

in  (33),  then  no  values  of  j  are  to  be  taken  in  (33). 

Since  [1,  Eq.  (12)]  is  used  for  p“,  the  factors  2,  4 j ,  2 j ,  and  4 

are  needed  in  [1,  Eqs.  (87) -(89)]  when  n  ^  1.  In  the  main  program,  lines 

48-51  store  these  factors  in  C3,  C4,  C5,  and  C6.  Line  53  stores  in  BB 

the  propagation  constant  of  the  medium  characterized  by  (y  ,  e  )  in 

[1,  Fig.  1],  This  propagation  constant  is  called  k  .  With  regard  to 

[1,  Eq.  (33)],  line  54  stores  nr  In  ET  and  line  55  stores  l/nr  in  ET1. 

Df  loop  23  puts  k  pj.  in  RA(I)  and  k  Zj  in  ZA(I) .  DO  loop  24  puts 

k  p  in  RB(I)  and  k  in  ZB(I). 

p  Q  i  P  Q  — 

The  elements  Y  and  Z  of  the  submatrices  on  the  right-hand 
nij  nij 

side  of  [1,  Eq.  (33)]  are  given  by  [1,  Eqs.  (36)  and  (35)].  The  testing 
pi  q  i 

functions  and  the  expansion  functions  ,TJ  appear  in  [1,  Eqs.  (36) 
and  (35)].  However,  the  submatrices  supplied  by  the  subroutine  YZA  are 

fg+  rs  t 

Y  and  Z  whose  elements  are  defined  by 
n  n 


42 


,.rs±  _ 

ni'j’  ~ 


zrs± 
ni' j ' 


h1<41'  •  °» 


r+  1  +  c  + 

o)> 
n 


t  ,4> 
t,4> 


(45) 


(46) 


Equations  (45)  and  (46)  are  (3)  and  (4)  with  the  superscript  ±  appended, 
with  i  replaced  by  i',  and  with  j  replaced  by  j ' .  In  (45)  and  (46),  the 

j*  +  si 

testing  functions  are  W  . ,  and  the  expansion  functions  are  J  . . .  To 

— ni  — nj 

pq  +  rs±  pq +  rs±  pi 

obtain  Y  ?.  from  Y'ti,  and  Z  ].  from  Z  , we  must  express  W  .  in 
nij  ni  j  nij  ni  j  m 

terms  of  and  J^J  in  terms  of  " , .  From  [1,  Eqs.  (23),  (24),  (27), 

(28),  (30),  and  (37)],  we  obtain 


wp:  =  (Jp7) 
— nl  -ni 


(47) 


where  *  denotes  complex  conjugate.  [1,  Eq.  (48))  and  (47)  state  that  the 

testing  functions  are  the  complex  conjugates  of  the  expansion  functions. 

pi  r± 

Hence,  the  expressions  for  in  terms  of  W^,  will  be  similar  to  the 

q+  gi 

expressions  for  2HJ  in  terms  of  , .  Consequently,  it  suffices  to  ex- 

q  i  si 

press  J  .  in  terms  of  J  ... 

-nj  -nj 

(1,  Eq.  (16)]  with  the  choice  of  superscript  +  is  rewritten  as 


1+M1A(2),  2+M1A ( 2 ) , . . .  L1A(2)  +  M1A(2) 


(48a) 

(48b) 

(48c) 

(48d) 


where 


MittiailMiHifa 
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J4*  + 

^n.j'-N  +2  ’ 


J'  -  1,  2,...  N  -2 


j'  =  N+  -  1,  N+,  ...  2N+-  3 


(49a) 


(49b) 


The  vector  functions  (49a)  are  the  t  directed  expansion  functions  on  the 
surface  (S+  +  A)  in  [1,  Fig.  4].  The  vector  functions  (49b)  are  the  <)> 
directed  expansion  functions  on  the  surface  (S+  +  A)  in  [1,  Fig.  4], 
Instead  of  expressing  each  j'  in  (48)  in  terms  of  j,  we  state  that 

~nj 

is  the  jth  vector  function  listed  on  the  right-hand  side  of  (48).  For 

example,  if  j  =  LlA(l)+3  and  if  L1A(2)>3,  then  J1+  =  J+  01W1A,0,. 

—  — nj  — n, 

In  (48), 


LlA(l)  =  Max(0,  M  -  1) 


L1A(2)  »  M 


L1A(3)  =  Max(0,  N+  -  M+  -  M  -  1) 


L1A(4)  =  N+  -  M+  -  M  -  LA 


MlA(l)  =  0 


M1A(2)  =  N  -2 


M1A(3)  =  M  +  M  -  1 


M1A(4)  =  N+  +  M+  +  M  -  3 


(50a) 


(50b) 


(50c) 


(50d) 


(51a) 


(51b) 


(51c) 


(51d) 


where  Max  denotes  maximum  value.  The  parameters  (L1A(I),  1=1, 2, 3, 4) 
and  (M1A(I),  1=1, 2, 3, 4)  are  calculated  by  lines  64-71  of  the  main 


program. 


The  vector  functions  (48a)  are  the  t  directed  expansion  functions 
on  the  first  part  of  the  outside  conductor.  The  first  part  of  the 
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outside  conductor  is  the  part  of  S+  below  the  aperture  in  [1,  Fig.  4].  The 
vector  functions  (48b)  are  the  <j>  directed  expansion  functions  on  the  first 
part  of  the  outside  conductor.  The  vector  functions  (48c)  are  the  t  directed 
expansion  functions  on  the  second  part  of  the  outside  conductor.  The  second 
part  of  the  outside  conductor  is  the  part  of  S+  above  the  aperture  in  [1, 

Fig.  4).  The  vector  functions  (48d)  are  the  $  directed  expansion  functions 
on  the  second  part  of  the  outside  conductor. 


J.,,  j'  -  1+M1B(1),  2+MlB(l), 

. . .  L1B(1)  +  M1B(1) 

(52a) 

JL,,  j’  =  1+M1B(2),  2+MlB(2) , 

. . .  LIB (2)  +  M1B(2) 

(52b) 

.T.,,  j'  =  1+M1B(3),  2+MlB(3), 

...  L1B(3)  +  M1B(3) 

(52c) 

JTj,,  j'  =  1+M1B(4) ,  2+M1 B ( 4 ) , 

. . .  L1B(4)  +  M1B(4) 

(52d) 

where 


,  j’  =  1,2,...  N  -2 


4~j'-N-+2  ’  J'  =  N__1’  N  2N"~3 


(53a) 


(53b) 


The  vector  functions  (53a)  are  the  t  directed  expansion  functions  on  the 
surface  (S  +  A)  in  [1,  Fig.  5].  The  vector  functions  (53b)  are  the  <t> 
directed  expansion  functions  on  the  surface  (S  +  A)  in  [1,  Fig.  5].  Equa¬ 
tion  (52)  means  that  J1,  is  the  jth  vector  function  listed  on  the  right- 


hand  side  of  (52).  In  (52), 

L1B(1)  -  Max(0,  M~  -  1)  (54a) 
L1B(2)  =  M"  ^  (54b) 
L1B(3)  -  Max(0,  N~  -  M_  -  M  -  1)  (54c) 
L1B(4)  -  N_  -  M“  -  M  -  LB  (54d) 
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M1B(1)  =  0  (55a) 
M1B(2)  =  tf  -  2  (55b) 
M1B(3)  =  M~  +  M  -  1  (55c) 
M1B(4)  =  N~  +  M_  +  M  -  3  (55d) 


The  parameters  (L1B(I),  I  =  1,2, 3,4)  and  (M1B(I),  I  =  1,2, 3, 4)  are  calcu¬ 
lated  by  lines  72-79  of  the  main  program. 


The  vector  functions  (52a)  are  the  t  directed  expansion  functions 
on  the  first  part  of  the  inside  conductor  and  (52b)  are  the  <J>  directed 
expansion  functions  there.  The  first  part  of  the  inside  conductor  is  the 
part  of  (S  +  A)  below  the  aperture  in  [1,  Fig.  5].  The  vector  functions 
(52c)  are  the  t  directed  expansion  functions  on  the  second  part  of  the 
inside  conductor  and  (52d)  are  the  4>  directed  expansion  functions  there. 
The  second  part  of  the  inside  conductor  is  the  part  of  (S  +  A)  above  the 
aperture  in  [1,  Fig.  5]. 

(1,  Eq.  (17)]  with  the  choice  of  superscript  +  is  rewritten  as 


,  j'  =  1+M2A(1),  2+M2A(l),  ...  L2(l)  +  M2A(1) 

,  j'  =  1+M2A(2) ,  2+M2A(2) ,  ...  L2(2)  +  M2A(2) 


(56a) 

(56b) 


+  2+ 
where  ,  is  given  by  (49).  Equation  (56)  means  that  is  the  jth 

vector  function  listed  on  the  right-hand  side  of  (56).  In  (56), 


L2(l) 


M+  t  0  ,  Ll A(4)  +  0 
M+  -  0,  Ll A(4)  +  0 
M+  i  0,  L1A(4)  =  0 
M+  *  0,  L1A(4)  =  0 


(57a) 
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L2(2)  =  M  -  1  -  LC 

(57b) 

M2A(1)  =  Max(0,  M+  -  1) 

(58a) 

M2A(2)  =  N+  -  2  +  M+ 

(58b) 

where  L1A(4)  is  given  by  (50d) .  The  parameters  L2(l),  L2(2),  M2A(1), 
and  M2A(2)  are  calculated  by  lines  80-85  of  the  main  program. 

The  vector  functions  (56a)  consist  of  the  t  directed  expansion 
function  which  straddles  the  first  part  of  the  outside  conductor  and  the 
aperture,  the  t  directed  expansion  functions  in  the  aperture,  and  the  t 
directed  expansion  function  which  straddles  the  aperture  and  the  second 
part  of  the  outside  conductor.  The  vector  functions  (56b)  are  the  (p 
directed  expansion  functions  in  the  aperture.  Four  different  cases  are 
necessary  in  (57a)  because  there  is  no  t  directed  expansion  function 
which  straddles  the  first  part  of  the  outside  conductor  and  the  aperture 
if  the  first  part  of  the  outside  conductor  is  absent  and  there  is  no  t 
directed  expansion  function  which  straddles  the  aperture  and  the  second 
part  of  the  outside  conductor  if  the  second  part  of  the  outside  conductor 
is  absent. 


[1,  Eq.  (17)]  with  the  choice  of  superscript  -  is  rewritten  as 

,  ,  j'  =  1+M2B(1),  2+M2B(l),  ...  L2(l)  +  M2B(1)  (59a) 

,  ,  J'  -  1+M2B(2),  2+M2B (2) ,  ...  L2(2)  +  M2B(2)  (59b) 

-  2- 
where  J  ,  is  given  by  (53).  Equation  (59)  means  that  J  is  the  jth 

~nj 

vector  function  listed  on  the  right-hand  side  of  (59).  In  (59),  L2(l) 
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and  L2(2)  are  given  by  (57).  M2B(1)  and  M2B(2)  are  given  by 

M2B(1)  =  Max (0,  M~  -  1)  (60a) 

M2B(2)  =  N~  -  2  +  M~  (60b) 

M2B(1)  and  M2B(2)  are  calculated  by  lines  86-87  of  the  main  program.  The 
vector  functions  (59a)  consist  of  the  t  directed  expansion  function  which 
straddles  the  first  part  of  the  inside  conductor  and  the  aperture,  the  t 
directed  expansion  functions  in  the  aperture,  and  the  t  directed  expansion 
function  which  straddles  the  aperture  and  the  second  part  of  the  inside 
conductor.  The  vector  functions  (59b)  are  the  <)>  directed  expansion  func- 

I 

tions  in  the  aperture.  Because  of  [1,  Eqs.  (15b)  and  (15d)],  the  vector 
functions  in  the  aperture  in  (59)  are  the  same  as  the  vector  functions  in 
the  aperture  in  (56). 

[1,  Eq.  (22)]  is  rewritten  as 
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calculated  by  lines  88-91  of  the  main  program.  The  vector  functions 
(61a)  are  the  t  directed  expansion  functions  in  the  aperture.  The  vector 
functions  (61b)  are  the  <f>  directed  expansion  functions  in  the  aperture. 

Replacement  of  the  superscript  +  in  [1,  Eq.  (22)]  by  the  super¬ 
script  -  gives 


v,j+M- 


,  j  =  1,2,  ...  M-2 


j37 

-nj 


-  h«  • J 


(64) 


Equation  (64)  is  consistent  with  [1,  Eqs.  (22),  (37),  (15b),  and  (15d)]. 
Equation  (64)  can  be  rewritten  as 


J  , ,  ,  j  '  =  1+M3B(1),  2+M3BC1),  ...  L3(l)  +  M3B(1) 
~tvj 


(65a) 


J37 

-nj 


J-.,  ,  j'  =  1+M3B(2),  2+M3 B ( 2 ) ,  ...  L3(2)  +  M3B(2)  (65b) 

“ nj 


3- 

where  J  is  given  by  (53).  Equation  (65)  means  that  J  .  is  the  jth 
— nj  —  nj 


vector  function  listed  on  the  right-hand  side  of  (65).  In  (65),  L3(l) 
and  L3(2)  are  given  by  (62a)  and  (62b).  M3B(1)  and  M3B(2)  are  given  by 

M3B(1)  =  M“  (66a) 

M3B(2)  =  N_  -  2  +  M_  (66b) 


M3B(1)  and  M3B(2)  are  calculated  by  lines  92-93  of  the  main  program.  The 
vector  functions  on  the  right-hand  side  of  (65)  are  the  same  as  those  on 
the  right-hand  side  of  (61) . 


To  accompany  the  expansion  functions  , 


of  (49)  and  (53),  testing 


4 

functions  W  ..  are  defined  by 
-ni 
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+  +  * 
Sti-  ‘  <£i*> 


(67) 


Equations  (49),  (53),  and  (67)  allow  (45)  and  (46)  to  be  recast  as 


<W' 


<W±.I 

ni 


,  ITQ*.  ,  0)> 

.  e^j;  ,  ,  o)> 

n  J 


(68) 

(69) 


Expressions  (68)  and  (69)  are  the  matrix  elements  supplied  by  the  sub¬ 
routine  YZA. 

Line  94  stores  in  JA  the  order  of  the  matrices  Y'+  and  Z+  whose 

n  n 

i'j'th  elements  are  given  by  (68)  and  (69).  Line  95  stores  in  JB  the 

order  of  the  matrices  Y'  and  Z  whose  i'j'th  elements  are  given  by 

n  n 

(68)  and  (69) .  Line  96  stores  in  IA  the  total  number  of  vector  func¬ 
tions  listed  on  the  right-hand  side  of  (48)  .  IA  is  the  order  of  the 
11+ 

submatrix  Z  in  [1,  Eq.  (33)].  Line  97  stores  in  IB  the  total  number 
n 

of  vector  functions  listed  on  the  right-hand  side  of  (52) .  IB  is  the 

order  of  the  submatrix  Z^1  in  [1,  Eq.  (33)].  Line  98  stores  in  N  the 

order  of  the  matrix  of  [1,  Eq.  (33)]. 

Line  99  puts  Y^+  of  (68)  in  YP  ( (n-Ml)*JA*JA+l)  to  YP ( (n-Ml+1) *JA*JA) 

and  Z+  of  (69)  in  the  corresponding  region  of  ZP  for  n  =  Ml,  Ml+1,  ...  M2. 

Storage  of  Y'+  and  Z+  in  YP  and  ZP  is  by  columns.  Line  100  puts  Y'  of 
n  n  v  n 

(68)  in  YM( (n-Ml) *JB*JB+1)  to  YM( (n-Ml+l)*JB*JB)  and  Z~  of  (69)  in  the 

n 

corresponding  region  of  ZM  for  n  =  Ml,  Ml+1,  ...  M2.  Storage  of  Y'  and 

n 

Z^  in  YM  and  ZM  is  by  columns.  With  reference  to  [1,  Eq.  (103)],  line  101 


puts 
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vte 

ni ' 

in 

R(i' 

+ 

Jl) 

i' 

=  1,2,. 

.  .NA-2 

V6 

ni' 

in 

R(i' 

+ 

J1  +  NA  -  2) 

i' 

=  1,2,. 

.  .  NA-1 

t<t> 

Ai' 

in 

R(  i ' 

+ 

Jl  +  JA) 

i’ 

=  1,2,. 

.  .NA-2 

V* 

ni 

in 

R(i ' 

+ 

Jl  +  JA  +  NA-2)  , 

i’ 

=  1,2,. 

. .NA-1 

(70) 


where 


J1  =  (n  -  Ml)  *  2  *  JA 


(71) 


and 


n  =  Ml,  Ml  +  1, 


M2 


(72) 


DO  loop  82  puts  k+Pj  in  RB(J)  .  k+Pj  is  needed  in  order  to  calculate 
[1,  Eq.  (88a)]. 

DO  loop  86  obtains  n  according  to  n  =  M-l.  If  n  ^  Ml,  DO  loop  89 

moves  Y'  +  down  into  YP(1)  to  YP(JA*JA)  and  Z+  down  into  ZP(1)  to 
n  n 

ZP(JA*JA),  DO  loop  90  moves  Y'“  down  into  YM(1)  to  YM(JB*JB)  and  z"  down 

n  n 

into  ZM(1)  to  ZM(JB*JB) ,  and  DO  loop  91  moves  the  V's  of  (70)  down  into 
R(l)  to  R(2*JA) . 

The  only  difference  between  Y^?T  of  [1,  Eq.  (36)1  and  Y'7,,.  of 

nij  ni  j 

(68)  lies  in  the  testing  functions  and  expansion  functions.  Likewise, 

nn  +  + 

the  only  difference  between  Z  . .  of  [1,  Eq.  (35)]  and  Z  of  (69) 

mj  m  j 

lies  in  the  testing  functions  and  expansion  functions.  The  expansion 
functions  in  [1,  Eqs.  (36)  and  (35)]  are  related  to  the  expansion  func¬ 
tions  in  (68)  and  (69)  by  (48),  (52),  (56),  (59),  (61),  and  (65).  Accord¬ 
ing  to  (67)  and  (47),  the  testing  functions  are  the  complex  conjugates 
of  the  expansion  functions.  Hence,  the  complex  conjugates  of  (48),  (52), 
(56),  (59),  (61),  and  (65)  relate  the  testing  functions  in  [1,  Eqs.  (36) 
and  (35)3  to  the  testing  functions  in  (68)  and  (69).  As  a  result. 
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Ypq± 

nij 


(73) 


P9±  =  71 
nij  ni' j ' 


Z 


(74) 


where  j'  is  the  subscript  of  the  jth  vector  function  listed  on  the  right- 
hand  side  of  either  (48),  (52),  (56),  (59),  (61),  or  (65),  whichever  is 
appropriate.  Similarly,  i'  is  the  subscript  of  the  ith  vector  function 
listed  on  the  right-hand  side  of  either  (48),  (52),  (56),  (59),  (61),  or 
(65),  whichever  is  appropriate.  For  example,  if  p  =  2  and  q  =  l  and  if  the 
superscript  +  is  chosen  in  (73),  then  (73)  becomes 


(75) 


where  j'  is  the  subscript  of  the  jth  vector  function  listed 
hand  side  of  (48)  and  i'  is  the  subscript  of  the  ith  vector 
on  the  right-hand  side  of  (56). 

With  regard  to  [1,  Eq.  (33)],  DO  loop  25  uses  (73)  and 

sn+ 

n 

0 

.21+ 


n 

31+ 


on  the  right- 
function  listed 

(74)  to  store 


(76) 


Ln  J 

by  columns  in  T.  The  elements  of  the  submatrices  in  (76)  depend  on  the 
expansion  functions  (48).  Hence,  j*  in  (73)  and  (74)  is  the  subscript  of 
the  jth  vector  function  listed  on  the  right-hand  side  of  (48).  The  index 
of  DO  loop  25  is  JJ.  If  JJ  «  1,  inner  DO  loop  26  obtains  the  values  of  j 
for  which  j'  is  given  by  (48a).  If  JJ  ■  2,  DO  loop  26  obtains  the  values 
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of  j  for  which  j*  is  given  by  (48b).  Similarly,  if  JJ  =  3,  DO  loop  26 

obtains  the  values  of  j  covered  by  (48c).  Finally,  if  JJ  =  4,  DO  loop 

26  obtains  the  values  of  j  covered  by  (48d). 

The  elements  of  Z^+  in  (76)  depend  on  the  testing  functions  W^. 

n  -ni 

These  testing  functions  are  related  to  ,  of  (67)  by  the  testing  func¬ 
tion  version  of  (48).  The  testing  function  version  of  (48)  is  (48) 
with  J^+  replaced  by  ,  J*  ,  replaced  by  vf^,,  and  j'  replaced  by  i' . 
It  is  now  apparent  that  i'  in  (74)  is  the  subscript  of  the  ith  vector 
function  listed  on  the  right-hand  side  of  the  testing  function  version 
of  (48).  Inside  nested  DO  loops  27  and  28,  line  141  stores  the  appro¬ 
priate  element  of  Z^1+  of  (76)  in  T.  The  index  of  DO  loop  27  is  II. 

If  II  =  1,  DO  loop  28  obtains  the  values  of  i  for  which  i'  is  given  by 
the  testing  function  version  of  (48a).  If  II  =  2,  DO  loop  28  obtains 
the  values  of  i  for  which  i'  is  given  by  the  testing  function  version 
of  (48b).  Similarly,  if  II  =  3,  DO  loop  28  obtains  the  values  of  i 
covered  by  the  testing  function  version  of  (48c).  Finally,  if  II  =4, 
DO  loop  28  obtains  the  values  of  i  covered  by  the  testing  function  ver¬ 
sion  of  (48d) . 

DO  loop  29  takes  care  of  the  0  in  (76). 

21+  2+ 
The  elements  of  Z  in  (76)  depend  on  the  testing  functions  W  . . 

n  “Til 

These  testing  functions  are  related  to  W^,  of  (67)  by  the  testing  func¬ 
tion  version  of  (56).  Hence,  i'  in  (74)  is  the  subscript  of  the  ith 
vector  function  listed  on  the  right-hand  side  of  the  testing  function 

version  of  (56).  Inside  nested  DO  loops  30  and  31,  line  157  stores 

21+ 

the  appropriate  element  of  Z  of  (76)  in  T.  The  index  of  DO  loop  30 

n 

is  II.  If  II  =  1,  DO  loop  31  obtains  the  values  of  i  for  which  i*  is 


i  - 
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given  by  the  testing  function  version  of  (56a).  If  II  -  2,  DO  loop  31 
obtains  the  values  of  i  for  which  i'  is  given  by  the  testing  function 
version  of  (56b). 

The  elements  of  Y"^+  in  (76)  depend  on  the  testing  functions  • 
n  — n  i 

These  testing  functions  are  related  to  ,  of  (67)  by  the  testing  func¬ 
tion  version  of  (61).  Hence  i'  in  (73)  is  the  subscript  of  the  ith 
vector  function  listed  on  the  right-hand  side  of  the  testing  function 

version  of  (61).  Inside  nested  DO  loops  32  and  33,  line  168  stores  the 

in¬ 
appropriate  element  of  Y^  of  (76)  in  T.  The  index  of  DO  loop  32  is  II. 

If  II  =  1,  DO  loop  33  obtains  the  values  of  i  for  which  i'  is  given  by 

the  testing  function  version  of  (61a).  If  II  ■  2,  DO  loop  33  obtains  the 

values  of  i  for  which  i'  is  given  by  the  testing  function  version  of  (61b). 

With  regard  to  [1,  Eq.  (33)],  DO  loop  34  uses  (73)  and  (74)  to  store 


n  zu~ 

r  n 

n  z21- 

r  n 
,,31- 


(77) 


by  columns  in  T.  Nested  DO  loops  34  and  35  obtain  the  values  of  j  in  (73) 

and  (74).  In  (73)  and  (74),  j'  is  the  subscript  of  the  jth  vector  function 

listed  on  the  right-hand  side  of  (52).  DO  loop  36  obtains  the  0  in  (77). 

Nested  DO  loops  37  and  38  obtain  the  values  of  i  for  of  (74)  in  which  i 

is  the  subscript  of  the  ith  vector  function  listed  on  the  right-hand  side  of 

the  testing  function  version  of  (52).  Nested  DO  loops  39  and  40  obtain  the 
21- 

values  of  i  for  Z  ..  of  (74)  in  which  i'  is  determined  by  the  testing  func- 
nij 

tion  version  of  (59).  Nested  DO  loops  41  and  42  obtain  the  values  of  i  for 
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Ynij  °f  in  which  i'  is  determined  by  the  testing  function  version 


of  (65). 


With  regard  to  [1,  Eq.  (33)],  DO  loop  43  uses  (73)  and  (74)  to 


store 


n  z12" 

r  n 

22+  n  _22- 
Z  +  d  Z 

n  r  n 

32+  32- 

Y  +  Y 

n  n 


by  columns  in  T.  Nested  DO  loops  43  and  44  obtain  the  values  of  j  in 

(73)  and  (74)  in  which  j'  is  determined  by  (56)  for  the  submatrices  with 

superscript  +  in  (78)  and  by  (59)  for  the  submatrices  with  superscript  - 

12+ 

in  (78).  Nested  DO  loops  45  and  46  obtain  the  values  of  i  for  Z  . .  of 

nij 

(74)  with  i'  determined  by  the  testing  function  version  of  (48).  Nested 

12- 

DO  loops  47  and  48  obtain  the  values  of  i  for  Z  . .  of  (74)  with  i'  de- 

nij 

termined  by  the  testing  function  version  of  (52).  Nested  DO  loops  49 

22+  22-  22+ 

and  50  obtain  the  values  of  i  for  Z  . .  +  n  Z  ,,  of  (78).  Z  ,  .  is  given 

nij  r  nij  nij  b 

22- 

by  (74)  with  i'  determined  by  the  testing  function  version  of  (56).  Z 

nij 

is  given  by  (74)  with  i'  determined  by  the  testing  function  version  of 

12+  i?_ 

(59).  Nested  DO  loops  51  and  52  obtain  the  values  of  i  for  Y  .  +  Y  “  of 

nij  nij 

32+ 

(78).  given  by  (73)  with  i’  determined  by  the  testing  function 

32- 

version  of  (61).  *s  given  by  (73)  with  i'  determined  by  the  testing 

function  version  of  (65). 

With  regard  to  [1,  Eq.  (33)],  DO  loop  53  uses  (73)  and  (74)  to 


store 
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_y23+  _  y23- 
n  n 

33+  1  33- 

+  (— )  ZJJ 
n  n  n 

r 


by  columns  in  T.  Nested  DO  loops  53  and  54  obtain  the  values  of  j  in 

(73)  and  (74)  with  j'  determined  by  (61)  for  the  submatrices  with 

superscript  +  in  (79)  and  by  (65)  for  the  submatrices  with  superscript  - 

13+ 

in  (79) .  Nested  DO  loops  55  and  56  obtain  the  values  of  i  for  °f 

(73)  with  i'  determined  by  the  testing  function  version  of  (48).  Nested 
DO  loops  57  and  58  obtain  the  values  of  i  for  Y^  of  (73)  with  i*  de¬ 
termined  by  the  testing  function  version  of  (52).  Nested  DO  loops  59 

9  9. 1.  in  n  i 

and  60  obtain  the  values  of  i  for  -Y  .  -  Y  .  of  (79).  Y  is  given 

nij  nij  nij 

by  (73)  with  i'  determined  by  the  testing  function  version  of  (56). 

23- 

Y  .  is  given  by  (73)  with  i'  determined  by  the  testing  function  version 
nij 

of  (59) .  Nested  DO  loops  61  and  62  obtain  the  values  of  i  for 

Z^+  +  (~)  of  (79).  is  by  (74)  with  i’  determined  by 

r  33- 

the  testing  function  version  of  (61).  Z  is  given  by  (74)  with  i' 

nij 

determined  by  the  testing  function  version  of  (65). 

Now  that  the  moment  matrix  T  has  been  stored  in  T,  we  turn  to  the 

n 

excitation  vector  B  of  [1,  Eq.  (101)].  The  V's  on  the  right-hand  side 
n 

of  [1,  Eq.  (103)]  reside  in  R.  Because  DO  loop  91  has  been  executed, 
storage  in  R  is  not  according  to  (70)  as  it  stands  but  according  to  (70) 


without  the  offset  Jl.  The  auxiliary  equation  involving  the  W's  in 
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[1,  Eq.  (103a)]  is  the  testing  function  version  of  (48).  Likewise,  the 
auxiliary  equation  in  [1,  Eq.  (103b)]  is  the  testing  function  version 
of  (56),  and  the  auxiliary  equation  in  [1,  Eq.  (103c)]  is  the  testing 
function  version  of  (61) . 

Based  on  the  considerations  in  the  preceding  paragraph,  lines 

0 

336-378  store  of  [1,  Eq.  (101)]  in  Y.  Nested  DO  loops  63  and  64 

i 4-A 

store  of  [1,  Eq.  (101)]  in  Y.  The  sign  factor  S  in  line  346 

compensates  for  the  minus  sign  in  front  of  in  (70).  The  number  of 

ni 

elements  in  the  row  vector  0  in  [1,  Eq.  (101)]  is  IB  because  this  row 

vector  can  be  traced  back  to  [1,  Eq.  (39b)].  DO  loop  65  stores  this 

row  vector  in  Y.  Nested  DO  loops  76  and  77  store  V  of  [1,  Eq.  (101)] 

n 

in  Y.  Line  364  compensates  for  the  minus  sign  in  front  of  V^,  in  (70). 

ni 

J  0 

Nested  DO  loops  66  and  67  store  1^  of  [1,  Eq.  (101)]  in  Y.  The  offset 

JA  in  line  370  is  mandated  by  the  superscript  0  in  in  [1,  Eq.  (103c)]. 

Line  376  accounts  for  the  net  effect  of  the  minus  sign  in  front  of  V^,  in 

[1,  Eq.  (103c)]  and  the  minus  sign  in  front  of  V^,  in  (70). 

Line  379  decomposes  T^  into  the  product  of  a  lower  triangular  matrix 

with  an  upper  triangular  matrix.  Line  380  stores  the  solution  X  to  [1, 

n 

Eq.  (65a)]  in  XX.  The  elements  [1,  Eq.  (68a)]  of  X^  are  the  I's  and  V’s 

n 

in  [1,  Eqs.  (87)-(89)].  Lines  384-385  store  e  in  Cl.  Line  386  stores 

n 

4j  in  C2.  If  L1A(2)  *  0,  then  it  is  evident  (50a)  and  (50b)  that  L1A(1)  =  0 
so  that  no  values  of  j  are  to  be  taken  in  (33)  and  (34).  If  L1A(2)  >  0, 
line  392  prints  (33)  and  (34).  If  L1A(4)  =  0,  then  it  is  evident  from 
(50c)  and  (50d)  that  L1A(3)  =  0  so  that  no  values  of  j  are  to  be  taken 
in  (35)  and  (36).  If  L1A(4)  >  0,  line  398  prints  (35)  and  (36).  Line 
403  prints  (37)  and  (38).  Line  409  prints  (39)  and  (40). 
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If  L2(l)  =  L2(2)  =  0,  then  it  is  evident  from  (57a),  (57b),  (62a), 

and  (62b)  that  L3(l)  =  L3(2)  =  0  so  that  no  values  of  j  are  to  be  taken 

in  (41)  -  (44).  If  either  L2(l)  >  0  or  L2(2)  >  0,  then  line  417  prints 

(41) — (42) .  If  L2(2)  =  0,  then  it  is  evident  from  (57b),  (62a)  and  (62b) 

that  L3(l)  =  L3(2)  =  0  so  that  no  values  of  j  are  to  be  taken  in  (43)  and 

(44).  If  1.2(2)  >  0,  line  421  puts  2j  in  Cl,  lines  422  and  423  put  2e  in 

n 

C2,  and  line  424  prints  (43)  and  (44). 


k 
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001  C  LISTING  OF  THE  MAIN  PROGRAM 

002  C  THE  SUBPROGRAMS  VZ A. BLOG .PLANE .DECOMP . SOLVE • 

003  //PGM  JOB  (XXXX.XXXX.2.2 >,' MAUTZ. JOE* ,REGI0N=275K 
0  04  //  EXEC  MATE*  V 
005  //GO. SYS  IN  OO  4 


AND  PRNT  ARE  NEEDED 


006  SJOB 
0  07 


008 
C09 
0  10 
Oil 


MAUTZ.T 1ME=5.PAGES=60 

COMPLEX  CI.C2.C3.C4.CS.C6. xx( 791. YPI2209)  . ZP C2209 I . Y M< 22091 
COMPLEX  ZM( 2209I.RC 2401 »T( 6241) .VI 79 > 

DIMENSION  XT  I I  0) .AT ( 10 ) . X(48 ) . A ( 48 ) .T  HR  13) .RA I 431 .ZAC  43)  .RB ( 43  > 
Ol  MENSION  28(431  «LIA(4).MI  A  C  4  )  .  LI  B  (  4  )  .  M  1  B  (  4  •  «L  2  (  2)  »M2A(  21  .  M.>H(  2) 
O IMENSION  L3( 2) .M3A(2J ,M3B( 2) . IPS( 79) 


012 

COMMON  C1.C2.XX.JX 

013 

REAO(l.lS)  NT.NPHI 

014 

IS 

FORMAT (213) 

015 

MR  I TE ( 3.9 )  NT.NPHI 

016 

9 

FORMAT ( *  NT  NPH l* / (X . I 3. I S ) 

0  1  7 

READ! 1 .10) 1 X  T 1 K ) . K= l • NT ) 

01  8 

READ! 1.I0)IAT<K).K=1 .NT) 

0  19 

10 

FORMAT 1 5E 14.7) 

020 

«RITE!3.I1)!XT!K).K=1 .NT) 

02  1 

MR  ITE(3.12)(AT(K) ,K=l.NT) 

022 

1  1 

FORMAT ( •  XT*/I 1X.5E14.7) ) 

323 

12 

FORMAT ( ■  Al* / ( 1X.5E14.7) ) 

024 

READ! 1. 10) ( X(K) .K=l .NPH I ) 

02S 

READ!  1.  10K  AIK)  .K=l.NPHl  ) 

026 

WRl  TE  (3 .13  >  (  X!  K)  .K  =  t  .  NPH  I  » 

027 

WRITE  (3.  141  (AIK)  ,K=1.NPHI  ) 

028 

13 

FORMAT! •  X*/( IX.SE14.7)) 

029 

14 

FORMAT ( •  A*/(1X.SE14.7)1 

03  G 

RE  AO ( 1.16)  NA.NB.MA.M8.MC.LA 

03  1 

16 

FORMAT! 10I3/4E 14.7) 

032 

MR ITE(3.17>  NA.NB. MA.M8.MC.L 

033 

1  7 

FORMAT! •  NA  NS  MA  MB  MC  LA 

034 

1  12X. *ER •  •  1  OX. • THR( 11 */!X.4E  l 

035 

READ (1.18)(RA( I ) . I— 1 . NA ) 

036 

READ! t.(S)(ZA( I ) . 1=1 .NA) 

037 

IB 

FORMAT! IOF0.4) 

038 

MRITEI3.19XRA!  (>•  I=1*NA) 

039 

MR! TE (3.20) IZAI I ) .1=1 .NA) 

040 

19 

FORMAT!  •  RAVI  IX.10F8.4)) 

041 

20 

FORMAT!*  ZA*2( IX. 10F8.4)) 

042 

RE AO! 1. IS) (RBI  I >. I=I.NB) 

043 

READ! 1 .18)128(1 1.1=1 .NB) 

C44 

MR  ITE!  3.21  MRS!  I  )  .1  =  1  .N8) 

045 

21 

FORMAT!*  Re*/I 1X.10FB.4)) 

046 

WR1TEI3.22  >!ZB(  0.1  =  1  .NB) 

04  7 

22 

FORMAT!'  ZB*^( 1 X.l OF0.4) ) 

048 

C3=2. 

049 

C4*4 . 4 1 0 • *  1 #) 

C50 

C 5=2.610. .1.) 

0S1 

C6=4. 

052 

UE=SQRT (UR4ER) 

053 

BB=BK*UE 

054 

ET=  SORT  (  UR  /ER  ) 

055 

ET 1=1 ./ET 

056 

OO  23  (=I.NA 

057 

R A! 1 )=BK«RA( 1) 

OSS 

ZA<! )*BK4ZA! 1 ) 

059 

23  CONTINUE 

060 

OO  24  1=1. NB 

LB  LC  Ml  M2*/IX.10I3/7X.*BK'.12X,*UR', 
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061 
06  2 
063 
064 
065 
066 
067 
068 
069 
070 
071 
072 
073 
074 
0  7S 
0  76 
077 
078 
0  79 
080 
081 
082 
083 
084 
085 
oeo 

097 
088 
089 
090 
091 
092 
093 
094 
095 
096 
097 
098 
099 
too 
I  01 
102 
103 
I  04 

105 

106 
IC7 
108 

109 

1 10 
I  1  I 
112 
I  13 
I  14 
I  IS 
I  16 
117 
I  18 

119 

120 


RBI  1 I =88* RSI  1 1 
ZB ( 1 >=00*Z0( II 
24  CONTINUE 

LI  All  )=MAXO(O.MA- l 1 
L1AI2I=MA 

LI  AOI^MAXu'O  «NAr— MA—MC— 1  I 

L  I  A<  4  )=NA-MA-MC— LA 

Ml  AC  I  1=0 

N1AI 2)=NA— 2 

NIA(3I>HA>MC-I 

Ml AI4  )=NA*MA*MC-3 

L181 1MMAX0I0.M8-1 I 

L 1B1 2 )=M8 

L 1 813 )=MAX0C  0* NB— MO— MC— 1 ) 

L 1 8 I 4 1 * NB-M0-MC-L8 
MI8I 1 1=0 
M18I2MN8-2 
Ml  813 )=MB*MC-I 
M1BI 4)>N8*MB+MC-3 
L2 1 1  l-MC 

IFIMA.EQ.O)  L2I1MMC-I 

IFCLlAI 4I.E0.0)  L2C 1I=L2( I  I  — l 

L  2 1 2  1*MC— 1— LC 

M2  A I  1  1=  MAX  01  0>  MA—  1  I 

M2 A ( 2 ) =NA— 24MA 

M28I 1 1=MAX0(0.MB-1) 

M28I2 )=NS— 2 AM8 
L3I1MMAX0 (0.MC-2I 
L3(  2)=L2(  21 
M3AI 1 )=MA 
M3A|2)=NA-2»MA 
M38I 1 1=M8 
M38I2 |=NB—24MB 
JA=2*NA— 3 
J8=2* N8— 3 

I  A=L  1  A(  1 1  4L IAC  2  )  *L  1  AC  31  *L  IA|  41 
IB=L1  81  1  ML  181 2  ML  181  3  ML  18141 
N=IA4-I8+L2{  1ML2C2ML3I  1ML3I2I 

CALL  YZA(M1*M2« NA < NPH I #NT •0»RA.ZA.X.A»XT«AT •  YP  «ZP) 
CALL  V2 A( Ml • M2  «NB  «NPH I »NT .O.R0.Z8.X. A. XT. AT t VM > ZMI 
CALL  PLANE I  Ml *M2 • I .NA.NT .RA, ZA» X T* AT .THR iRI 
DO  82  J= I . N8 
RBI J)=R8I JI/UE 
82  CONTINUE 
J  A2=  J  A*  J  A 
J 82=38* J8 
M3*MM1 
M4=M2*| 

OO  86  M=M3*M4 
IFfM.EU.M3)  GO  TO  88 
J2=M-M3 
•ft  =  J24  JA2 
OO  89  J=I.JA2 
VP(3)=VP(3» Jll 
ZPI Jl-ZPI J*J1 I 
89  CONTINUE 
JI*JI»J02 
OO  90  3= 1 .  J82 
YMIJI«VM(3»Jt I 
ZMI Jt-ZMI J*JI  | 
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121 

90 

CONTINUE 

122 

J3=2*JA 

1  23 

Jl=J2* J3 

124 

DO  91  J*I.J3 

125 

R1 J»  =  RI J4J11 

126 

91 

CONTINUE 

127 

88 

JT»0 

1  20 

OO  23  JJ*1,4 

129 

J2=L1A< JJ1 

130 

IF< J2-EQ.0)  GO 

TO 

25 

131 

J1=MIAI JJI4JA 

1  32 

OO  26  J*1  .  J2 

133 

OO  27  0*1.4 

134 

1 3=L1 A(  X  I  ) 

1  35 

IF! 13.60.0)  GO 

TO 

27 

1  36 

J3=J1  +  M  IA< 1 I ) 

137 

I l=J34l 

133 

I2-J3fl3 

139 

DO  28  1*11*12 

1  40 

JT=JT4l 

141 

T ( JT1*ZPC 1 1 

142 

28 

CONTINUE 

143 

27 

CONT (NUE 

I  44 

IFOB.EQ.O)  GO 

TO 

79 

1  45 

OO  29  1*1.18 

1  46 

JT=JT41 

147 

T<  JTl-O. 

1  48 

29 

CONTINUE 

1  49 

79 

OO  30  I 1*1 *2 

1  50 

13*1.2(111 

151 

1F(I3.EO.O)  GO 

TO 

30 

152 

J3=JI  4M2AC  11  ) 

1  53 

I  1=J34I 

154 

I  2=  J3H3 

155 

OO  31  1*11*12 

156 

JT*JT41 

157 

T( JT»*ZPCI) 

158 

31 

CONTINUE 

1  59 

30 

CONTINUE 

l  60 

DO  32  I 1*1 *2 

161 

13*1.301) 

162 

IF(U.EO.O)  GO 

TO 

32 

163 

J3*Jl 4M3A( II) 

1  64 

1 l*J3*l 

165 

12= J3+ I 3 

166 

OO  33  1*11 *12 

1  67 

JT*JT4I 

168 

T(JT)*YPI 1 ) 

169 

33 

CONTINUE 

1  /  C 

32 

CONTINUE 

1  71 

JI*JI»JA 

1  72 

26 

CONTINUE 

l  73 

25 

CONTINUE 

1  74 

00  34  J J* 1*4 

1  75 

J2*LI8( J J) 

1  76 

IF ( J2*E0*  0 1  GO 

TO 

34 

1/7 

J1*M|B( J J ) ♦ J8 

1  78 

00  35  J*I*J2 

l  79 

IFIIA.EQ.OI  GO 

TO 

80 

1  80 

OO  36  I" I*  I A 

\ 

tat 

182  T(JI)*0. 

183  36  CONTINUE 

184  80  OO  37  11=1.4 

185  1 3=L  1  U(  II) 

186  IF1I3.EQ.0)  CO  TO  37 

187  33=3|4MlB(lt) 

188  I|sj3it 

139  12=334(3 

190  OO  38  1=11.12 

191  3T=3T4l 

192  T ( JT  >=£T  *2N( 1 ) 

193  38  CONTINUE 

194  37  CONTINUE 

195  OO  39  1 1=1 .2 

196  1 3=1-2  III) 

197  IF! (3.EQ.0)  GO  TO  39 

198  33=314M2B< Ill 

199  11=3341 

200  12=33413 

201  OO  40  1=11.12 

202  JT=JT*1 

203  T< JT»=ET42M1 I » 

204  40  CONTINUE 

2  OS  39  CONTINUE 

206  OO  41  11=1.2 

20 7  I3=L3T 111 

208  IF<  O.EO.O)  GO  TO  41 

209  33=3 I 4N3B1 111 

210  11=3341 

211  12=33413 

212  OO  42  1=11.12 

213  3T=3T41 

214  M3T)=YMI1) 

215  42  CONTINUE 

216  41  CONTINUE 

217  3I=J|43B 

218  35  CONTINUE 

219  34  CONTINUE 

220  OO  43  33=1*2 

221  J  2=1-2  (33) 

222  (F132.EQ.0)  GO  TO  43 

22 3  31A=M2A< 33)434 

224  3IB=M2B(J3)*3S 

225  DO  44  3=1 *32 

226  OO  45  11=1.4 

227  I 3=L 1  A<  1 1 1 

228  IF113.EO.O)  GO  TC  45 

229  33=31 A4WI A 1 1 1 ) 

230  11=3341 

231  12=334(3 

232  OO  46  (=11.12 

2J3  3T=3T  41 

234  TlJT)-ZP(l) 

235  46  CONTINUE 

236  45  CONTINUE 

237  OO  47  (1=1*4 

238  I  3=1.  IBl  1 1  ) 

239  IFI 13*E0*0)  GO  TO  47 

240  33=31  S4M1  Bill) 


62 


241  Il=334| 

242  12=33413 

24  J  DO  48  1=11.12 

244  3T=3T4l 

245  T  (  3T  J=£  T4ZMI I ) 

246  48  CONTINUE 

247  47  CONTINUE 

248  OO  49  I  1  =  1  .2 

249  13=L2(II1 

250  IFII3.E0.0)  GO  TO  49 

251  J3  =  J1  A4M2AOI) 

252  I 1=3341 

253  12=33413 

254  J4=JIB4M2B< I  II 

255  OO  SO  (=11.12 

256  J4 - 34 ♦ 1 

207  3T=3T4l 

258  T( JT»=ZP< I )4ET«ZM< J4» 

259  50  CONTINUE 

260  49  CONTINUE 

261  OO  St  l  (= I *2 

262  I3=L3 (11) 

263  1FI13.EQ.0)  GO  TO  51 

264  33=J1A4M3A< II ) 

265  1 1=331-1 

266  I2=J34(3 

267  J  4= J 1 84M30 (II) 

268  OO  S2  1=11.12 

269  34=3441 

2  70  J  T=3T4 1 

271  TIJT)=YP(I 14YMI34) 

272  52  CONTINUE 

273  51  CONTINUE 

2  74  3  1  A=  J  1  A43  A 

275  31E=J1B43B 

2  76  *4  CONTINUE 

277  43  CONTINUE 

278  OO  53  33=1.2 

279  32=L3<33> 

280  I F 1 32 . £0.0 )  GO  TO  53 

221  3IA=M3A< 33>*3A 

282  31B=M3B(33 )*3B 

283  OO  54  3=1.32 

284  DO  55  11=1.4 

285  13=1.1  AII11 

286  l Ft  13. EQ. 0 }  GO  TO  55 

287  33=31A4MIA(I1> 

288  11=3341 

289  12=33413 

290  OO  56  1=11.12 

291  3  T=  J  T  4  1 

292  T<JT»=-VPin 

293  56  CONTINUE 

294  SS  CONTINUE 

295  OO  57  11=1.4 

296  1 3=L 1 B1 1 1 1 

297  IFII3.EO.O)  GO  TO  57 

298  33=  J104M18IIIJ 

299  11=3341 

300  12=33413 


63 


301 

302 

303 

304 

305 

306 

307 
306 

309 

310 
3  1  1 

312 

313 

314 

315 

316 
3  1  7 
3  I  8 
3  19 
320 
32  1 

322 

323 

324 

325 

326 

327 

328 

329 

330 

331 

332 

333 

334 

335 
3  36 

337 

338 

339 

340 

341 

342 

343 

344 

345 

346 

347 

348 

349 

350 

351 

352 

353 

354 

355 

356 

357 

358 

359 

360 


90  58  1=11.12 

JT=3T4t 

T<  3T)=-YM(  1) 

58  CONTINUE 
57  CONTINUE 

OO  59  (  1=1  .2 
1 3=1.2  <  I  1) 

IF! 13.EQ.0)  GO  TO  S9 

J3= Jl A*  M2 At  1(1 

11  =  3341 

12=33413 

J4=J 1 8»M28 ( 111 

OO  60  1=11 .12 

34=3441 

3T=3T  4 1 

T< JT1=-VP( I1-YM1341 

60  CONTINUE 

59  CONTINUE 

OO  61  11=1.2 

13=1.31  I  II 

IF113.EQ.01  GO  TO  61 

3  3=3 1 A  4M3A 1 1 1 1 

11=3341 

12=33413 

34=3184  8381  111 

OO  62  1=1 1,12 

34=34  41 

3T=3T41 

T(3TJ=ZPC  I  )4ET1«ZM<  341 

62  CONTINUE 

61  CONTINUE 
31A=31A43A 
318=318438 

54  CONTINUE 
53  CONTINUE 
3Y=0 
S=l« 

OO  63  11=1.4 
13=1.1  A1  11) 

IF113.EQ.01  GO  TO  83 
33=M|A1I1) 

1 1=3341 

12=33413 

OO  64  (=11 .12 

3V=3V41 

V 1 3 Y 1=  S4R ( I ) 

64  CONTINUE 
83  S=-S 

63  CONTINUE 
IF118.EQ.0)  GO  TO  81 
OO  65  (=1.18 

3Y=3Y  41 
Y  1  3Y  1=0  . 

65  CONTINUE 

81  OO  76  11=1.2 
13=1.2(1  I) 

(F ( 13. EQ. 0)  GO  TO  76 
33=M2A( 1 1 ) 

11=3341 

12=334(3 


23S 


2 
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36  K 

OQ  77  1=1 1.(2 

362 

JY=  JYM 

3*3 

VI JY)=R(I) 

364 

IF(!(.EQ.2>  Y(JY)=-Y<JY) 

365 

77 

CONTINUE 

366 

76 

CONTINUE 

36  7 

OO  66  I  1  =  1  .2 

360 

I  3=L3(  11) 

369 

IF  1  I  3.EQ.0 >  GO  TO  66 

3/0 

J3=M3A< I I )  ♦  JA 

371 

1 1= JJ» 1 

372 

I2=J34I3 

37  3 

OO  67  1  =  11 .12 

3  7  4 

JY=JTM 

3/5 

Y1 JY)=R1 I) 

3  /6 

IF1II.EQ.2)  Y1JYI=-Y1JY) 

3/7 

67 

CONTI NUE 

378 

66 

CONTINUE 

379 

CALL  OECOMPlN. IPS.T ) 

380 

CALL  SOLVE  IN. IPS.T.Y.XX I 

38  1 

MS=M-  1 

382 

MRITEI3.96)  MS 

383 

96 

FORMAK'O*  .  13.*  TH  MODE  ELECTRIC  ANO  MAGNETIC  CURRENTS*  > 

384 

C1=C3 

385 

IF1M5.EQ.0)  r 1=1. 

386 

C2=C4 

387 

JX=0 

388 

J2=H  AI2) 

389 

IF1J2.E0.0)  GO  TO  84 

390 

MRI TE13.68) 

391 

68 

F 0  RM  AT ( • OEL  ECTR IC  CURRENT  ON 

FIRST  PART  OF 

OUTSIDE 

CONDUCTOR* ) 

392 

CALL  PRNT1LI A< 1 1.J2.1. l.RA) 

393 

84 

J2=L 1 Ai 4) 

394 

IF1J2.EQ.0)  GO  TO  74 

295 

WRITE13.69) 

396 

69 

FORMAT! *OELECTRIC  CURRENT  ON 

SECOND  PART  OF 

OUTS  IDE 

CONDUCTOR*) 

297 

J3=MA4MC 

298 

CALL  PRNT(L1A(3). J2,  J3.  J3.RA) 

399 

74 

J2=H  B12) 

400 

I F ( J2.EO.O )  GO  TO  85 

40  1 

WRITE (3.70) 

402 

70 

FORMAT1»OELECTRIC  CURRENT  ON 

FIRST  PART  OF 

INSIDE  CONDUCTOR*) 

403 

CALL  PRNTILI81 1) .J2.1.1.RS) 

404 

85 

J2=LIG(4) 

405 

IF1J2.E0.0)  GO  TO  75 

406 

MRI TE 1 3. 71) 

407 

71 

FORMAT I'OELECTRIC  CURRENT  ON 

SECOND  PART  OF 

INSIDE 

CONDUCTOR* ) 

408 

J3=MSfrMC 

409 

CALL  PRNT (L1B(3).J2.J3.J3.R8) 

4  10 

75 

JI=L2(I  1 

41  1 

J2=L2(2) 

4  12 

I F 1 1 31FJ2) •£ Q. 0)  GO  TO  86 

4  13 

MR  1 TE (3*72 ) 

41  4 

72 

FORMAT! 'OELECTRIC  CURRENT  IN 

APERTURE* ) 

4  15 

J J=M2A( 1)41 

416 

J4=MAM 

41  7 

CALL  PRNTI Jl. J2.J3. J4.RA) 

418 

IF1J2.E0.0)  GO  TO  86 

4  1  9 

MRITE13. 73) 

420 

73 

FORMAT!* OMAGNET IC  CURRENT  IN 

APERTURE* ) 

r  .  . _  yj 


421  C1=C5 

422  C2^C6 

423  IFIMS.EQ.OI  C2=2. 

424  CALL  PRNT (L3(  1  )  . J2. J4, J4. RA ) 

425  86  CONTINUE 

426  STOP 

427  ENO 

*OATA 

2  20 

-0 .57735036700  0.5773503E700 
0*1 OOOOOOE  V  01  0. 1 ooooooe*o 1 

-0.9931286E700-0.9639719E  70 0- 0. 9 1 22 344E700- 0 . 839 1  l 706700-0 .74633 l 9E 700 
-0 .6360 53 7E *00-0. 51 086 70E 700-0. 3 7 3706 l E *0 0- 0 .2 2 778596 7 00 -0 . 76S2652E-0 I 
O. 7652652E-01  0.22778S9E700  0 .3 73 706 1 E *00  O .5  I  086 70 E 700  0  .636053 7E  +  00 

0.746331 9E  7  00  0.83911706700  0 . 9  I  223 44E 7 0 0  0. 963 97 1 9E700  0 .993  I 286E 700 
0.1 761401E-0I  0.4060 143E-01  0. 6  26 72 05E- 0 1  0. 832 76 74E- 0 1  0 . 1 0 1 930  IE 700 
0. 11 81 945E700  O. 13 1 6886E700  0. 1420961E700  0 . 1 4 9 1 7 30 E 700  0 . 1 527534E 700 
0. 1527S34E 700  0. 1 49 1 730E 7 00  0.  1 4 2096 1 E7 00  0. 1 31 6886E700  0 . 1  1 8 1945E 700 
0 . 10 19301E700  0.8J27674E-01  0 . 626 7205E- 0 1  0. 4 06 0 1 43E-01  0. 1 761 4 01 E-0 1 

1513  54  5  000  1  1 

0.2500000E70 1  0.1000000E701  0 . 2000000E7 0 1  0. OOOOOOOE 700 


0. 0000 

0.2588 

0.5000 

0.7071 

0.8660 

0.9659 

0.9914 

t .0000 

0.9914 

0.9689 

0.  8660 

0.7071 

0. £000 

0.2588 

0.0000 

-1.0000 

-0.9659 

-0.8660 

-0.7071 

-0.5000 

-o.2sea 

-0.  1305 

0.0000 

0.1305 

0.2588 

0.5000 

0.7071 

0.8660 

0.9659 

1.0000 

0.0000 

0.2415 

0. 4830 

0.7244 

0. 9659 

0.9914 

1.0000 

0.9914 

0 .9659 

0.7244 

0.4830 

0.24  IS 

0.0000 

-0.2588 

-0.2588 

-0.2588 

-0.2588 

-0.2588 

-0.1305 

0.0000 

0.1305 

0.2588 

0.2388 

0.2588  0.2568  0.2588 

SSTOP 
/* 

// 

PRINTED  OUTPUT 
NT  NPHI 
2  20 
XT 

-0.5773S03E700  0.5773503E 700 
AT 

0. 1000000E70I  0. 1000000E701 
X 

-0.9931 2 86E 7  00-0. 96  3  97 19E700-0 . 9 12234 4E 700-0. 8391  I 70E7 00-0 . 74633 19E 70  0 
—  0 .6360537E  7  00—0.51 08670E  7  00—  0. 3  737  06 1 E7  00— 0 . 22  77859E700— 0 .7652652E-01 
0 • 7652652E— 0 1  O .227 7859E 700  0 . 373 706 l £7 0 0  0. 5 1 08670E700  0. 6360S3 7E700 
0. 746331 9E700  0.8391 1 70E700  0 .9 1 223446700  0.96397 19E700  0 .9931 286E700 

A 

0.1 76I401E-0 1  0.4060 I43E-01  0. 626 72 08E-0 I  0. 83 27675E- 0 l  0. I  01 9301 E 700 
0.1 181945E700  0.131 6 8 86 E 700  0 . 1 42096  I E 700  0. I  4 9 1 730E700  0 . 1 52 T534E 700 
0.15275346  700  O. 1 49 1 730E 7 00  0.  1 42096 IE  700  O. 1 31 6886E700  0 . 1 1 8 1 945E 700 
0 . 10193016700  0.83276 7SC— 01  0 .626 72 0 8E- 0 1  0.4060 1 43E-01  0. 1 76 1 40 1E-0 1 

NA  N8  MA  MB  MC  LA  LB  LC  Ml  M2 
15  13  545000  1  1 

DK  UR  ER  T  HR {  1 ) 

0.2  5000  006  7  01  0. 1 000000  £70  I  0 . 2000000E 701  0 . 0 OOOOOOE 700 

RA 

0.0000  0.2588  0.5000  0.7071  0.0660  0.9659  0.9914  1.0000  0.9914  0.9659 

0.8660  0.7071  0.5000  0.2588  0.0000 

ZA 


UP 


0.  0000  0.2  V I  5  0.4830  0.7244  0.9659  0.9914  1.0000  0.9914  0.96S9  0.7244 

0.4810  0.2415  0.0000 

Ztl 

-0.2588  -0.2*508  -0.2588  -0.2588  -0.2588  -0.1305  0.0000  0.1305  0.2588  0.2588 

0.2588  0.2533  0.2588 

1  TH  MODE  ELECTRIC  A NO  MAGNETIC  CURRENTS 
ELECTRIC  CURRENT  ON  FIRST  PART  OF  OUTSIDE  CONDUCTOR 

REAL  Jf  (MAG  JT  MAG  JT 

0.  1. 1?  ST.  t  00-0.3  7  l  5  E  1 0  0  0.1077E401 
O  .  3889E  tOO- 0.322 J£- 01  0.5946E* 00 
0.2  590E  tOO  0.7627E»00  O.GOSSEtOO 

—  0.  42  OSE  t  00  0.1188Et01  0 . 1 26 1 E tO 1 

RFAL  JP  IMAG  JP  MAG  JP 

-0.60926t00  O. 1 065E  1 0 1  0.1227Et01 

—  O  •  5  .1  8  6E  tOO  0.6909EE00  0.8760Et00 
-0.3J-*8EtG0  0.41  60E  t  00  0.5J39Et00 

—  0*80  79E  —  01  0.2  986  E  1 00  0.3074£t00 

-0.S854E-01  0.530  8E  *00  0.5340£t00 

ELECTRIC  CURRENT  ON  SECOND  PART  OF  OUTSIDE  CONDUCTOR 

REAL  JT  IMAG  JT  MAG  JT 

-O  .«426Et0O-0-l79IE»0  1  0.1844Et01 
O.  4  3  00E«-00-0.  I  9  04  £  to  1  0.1957Et0l 
0.1l05Et01-0.lft07E*Ol  0- I950£t0t 
0.!478Et01-O.1261Et01  0-194JEt0l 

REAL  JP  IMAG  JP  MAG  JP 

-0. 394  IE t 00- 0.53 44 Et 00  0.1042Et0l 
-0.2l90Et0O-0.127IEt01  0.1290EF01 
0  .5212£t00-0. I483£t01  0.1572E  +  01 

0. 1 l l 3ET01— 0.l3ft5Et0I  0.1765Et01 
0. I525£t0l-0. 1 I836t0l  0.1930Et01 

ELECTRIC  CURRENT  ON  FIRST  PART  OF  INSIDE  CONDUCTOR 

REAL  JT  IMAG  JT  MAG  JT 

O.I267E»OI  0.38  71 E— 01  O.1268Et0I 
0 • 26  08E  tOO  0.5877E-0I  0.2673EF00 
-0.7252E400  0.1 498E tOO  0.7415E400 

REAL  JP  IMAG  JP  MAG  JP 

-0.1 558E to 1  —  0 .3328E—  01  0. 1558E tOl 
-0.11166*01-0 .6796E-01  0.1 1 ISEtO t 

-O.70<S2Et0O-O.762IE-0I  0.71  03Et00 
-0.1 70'Et 00-0. 31 73Et CO  0. 3603E tOO 

ELECTRIC  CURRENT  ON  SECONO  PART  OF  INSIOE  CONDUCTOR 

REAL  JT  IMAG  JT  MAG  JT 

-0. 75 49E tOO— 0.2851 EtOO  0.8005EtOO 
0.251  3F.  tOO—  0.485  7E—  0 1  0.2S56E«-00 
0.1200E+0I  0.8836E-0I  0-1283Et0l 

REAL  JP  IMAG  JP  MAG  JP 

0.657  7E  tOO— 0 .I862£t00  0.6932E+00 
0 • 7V89E  tOO  0 . 742  BE—  02  0.7989E»00 


\ 

C- 


0.1129E*01  0.67S1E— 0 1  0.1l31E*0l 

0.1S84E*01  0*121 4E  *00  0.15896*01 

ELECTRIC  CURRENT  IN  APERTURE 

REAL  JT  I  MAG  JT  MAC  JT 

-0.11026*01  0  *6  /  706+  00  0. 12946*01 
—0 • 12616*01  0.3286E*00  0.13036*01 
-0. 1 339E*  01  —  0. 4550E— 0 1  0-1 339  6*01 

—  0  »  I  3 056 *0 1— O •  44  786*  00  0.  13806*01 
-0.I196£*01-0 .9 0236*00  0.14986*01 

REAL  JP  l  MAG  JP  MAG  JP 

— 0 • 1 580E*00  0.26586*00  0. 30926*00 

—  0-2  7286*00  0* I  5346*00  0.31306*00 

—  O .3 7S0E *00  0*  5005E—  0 1  0. 3783E*00 
-0 .5571 E *00— 0 »5482E— 0 1  0.55986*00 

MAGNETIC  CURRENT  IN  APERTURE 

REAL  JT  1 MAG  JT  MAG  JT 

0.14296*00  0.7S46E-01  0.16166*00 
0.15896*00  0.13936*00  0.21136*00 
0.10496*00  0. 1 5 79E*  00  0.  1896E*00 

REAL  JP  I  MAG  JP  MAG  JP 

0.65206*00  0.154 16*00  0.6708E*00 
0.1SS9E*00  0.38556*00  0.41S8E*00 
-0.10006*00  0.45236*00  0.46336*00 
-0.83956*00  0.48116*00  0.96756*00 
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R  ionospheric  propagation,  solid  state  sciences,  microwave  ** 

physics  and  electronic  reliability,  maintainability  and 
s  compatibility , 
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